# HG changeset patch # User John Donoghue # Date 1578309576 18000 # Node ID aad0cc165a27bc9e96db41191e6ca9fc55b21506 # Parent 74fe6c154dee07355b2a2f91afc02c01ae583def of-communications: update to v1.2.2 * dist-files.mk: remove refs to patches * src/of-communications.mk: update version, checksum * src/of-communications-1-fixes.patch, src/of-communications-2-fixes.patch, src/of-communications-3-fixes.patch, src/of-communications-4-fixes.patch, src/of-communications-5-fixes.patch, src/of-communications-6-deprecated.patch : removed files diff -r 74fe6c154dee -r aad0cc165a27 dist-files.mk --- a/dist-files.mk Wed Jan 01 18:52:37 2020 -0500 +++ b/dist-files.mk Mon Jan 06 06:19:36 2020 -0500 @@ -451,12 +451,6 @@ ocaml-native-1-fixes.patch \ ocaml-native.mk \ of-actuarial.mk \ - of-communications-1-fixes.patch \ - of-communications-2-fixes.patch \ - of-communications-3-fixes.patch \ - of-communications-4-fixes.patch \ - of-communications-5-fixes.patch \ - of-communications-6-deprecated.patch \ of-communications.mk \ of-control.mk \ of-data-smoothing.mk \ diff -r 74fe6c154dee -r aad0cc165a27 src/of-communications-1-fixes.patch --- a/src/of-communications-1-fixes.patch Wed Jan 01 18:52:37 2020 -0500 +++ /dev/null Thu Jan 01 00:00:00 1970 +0000 @@ -1,78 +0,0 @@ -# HG changeset patch -# User John Donoghue -# Date 1441300743 14400 -# Thu Sep 03 13:19:03 2015 -0400 -# Node ID eeb20b3c6074e8932e40b14c9f99985340ceeca5 -# Parent dedf92d61d33525da3f71a635979c447defb438a -Check for ls-oct-ascii.h, ls-oct-text.h (Bug #45856) - -src/configure.ac: check for ls-oct-ascii.h, ls-oct-text.h and set HAVE_XXXX variables accordingly. - -src/ov-galois.cc: include either ls-oct-ascii.h, ls-oct-text.h dependong on HAVE_XXXX variables - -diff -r dedf92d61d33 -r eeb20b3c6074 src/configure.ac ---- a/src/configure.ac Sat Apr 04 12:20:33 2015 -0400 -+++ b/src/configure.ac Thu Sep 03 13:19:03 2015 -0400 -@@ -99,5 +99,21 @@ - PKG_CPPFLAGS="$PKG_CPPFLAGS -DHAVE_OCTAVE_BASE_VALUE_PRINT_CONST=1" - fi - -+comm_save_CXXFLAGS=$CXXFLAGS -+CXXFLAGS="$CXXFLAGS $comm_CXXFLAGS" -+AC_LANG_PUSH(C++) -+AC_CHECK_HEADER([octave/ls-oct-ascii.h], [PKG_CPPFLAGS="$PKG_CPPFLAGS -DHAVE_OCTAVE_ASCII_H=1"], [], -+[ -+#include -+#incude -+]) -+AC_CHECK_HEADER([octave/ls-oct-text.h], [PKG_CPPFLAGS="$PKG_CPPFLAGS -DHAVE_OCTAVE_TEXT_H=1"], [], -+[ -+#include -+#include -+]) -+AC_LANG_POP(C++) -+CXXFLAGS=$comm_save_CXXFLAGS -+ - AC_CONFIG_FILES([Makefile]) - AC_OUTPUT -diff -r dedf92d61d33 -r eeb20b3c6074 src/ov-galois.cc ---- a/src/ov-galois.cc Sat Apr 04 12:20:33 2015 -0400 -+++ b/src/ov-galois.cc Thu Sep 03 13:19:03 2015 -0400 -@@ -30,7 +30,12 @@ - #include - #include - --#include -+ -+#if defined(HAVE_OCTAVE_TEXT_H) -+# include -+#else -+# include -+#endif - - #include "galois.h" - #include "ov-galois.h" -# HG changeset patch -# User John Donoghue -# Date 1448402937 18000 -# Tue Nov 24 17:08:57 2015 -0500 -# Node ID e8b6bae07002b9e94bfc104b1e844963ccaa5f41 -# Parent eeb20b3c6074e8932e40b14c9f99985340ceeca5 -define OV_REP_TYPE it it is not defined (Bug #46521) - -* src/ov-galois.h: define OV_REP_TYPE if not defined already - -diff -r eeb20b3c6074 -r e8b6bae07002 src/ov-galois.h ---- a/src/ov-galois.h Thu Sep 03 13:19:03 2015 -0400 -+++ b/src/ov-galois.h Tue Nov 24 17:08:57 2015 -0500 -@@ -44,6 +44,10 @@ - #endif - #endif - -+#if ! defined (OV_REP_TYPE) -+# define OV_REP_TYPE octave_base_value -+#endif -+ - class octave_value_list; - class tree_walker; - diff -r 74fe6c154dee -r aad0cc165a27 src/of-communications-2-fixes.patch --- a/src/of-communications-2-fixes.patch Wed Jan 01 18:52:37 2020 -0500 +++ /dev/null Thu Jan 01 00:00:00 1970 +0000 @@ -1,399 +0,0 @@ -diff --git a/src/base-lu.cc b/src/base-lu.cc -new file mode 100644 ---- /dev/null -+++ b/src/base-lu.cc -@@ -0,0 +1,191 @@ -+/* -+ -+Copyright (C) 1996-2015 John W. Eaton -+Copyright (C) 2009 VZLU Prague -+ -+This file is part of Octave. -+ -+Octave 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 3 of the License, or (at your -+option) any later version. -+ -+Octave 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 Octave; see the file COPYING. If not, see -+. -+ -+*/ -+ -+#ifdef HAVE_CONFIG_H -+#include -+#endif -+ -+#include "base-lu.h" -+ -+template -+base_lu::base_lu (const lu_type& l, const lu_type& u, -+ const PermMatrix& p) -+ : a_fact (u), l_fact (l), ipvt (p.transpose ().col_perm_vec ()) -+{ -+ if (l.columns () != u.rows ()) -+ (*current_liboctave_error_handler) ("lu: dimension mismatch"); -+} -+ -+template -+bool -+base_lu :: packed (void) const -+{ -+ return l_fact.dims () == dim_vector (); -+} -+ -+template -+void -+base_lu :: unpack (void) -+{ -+ if (packed ()) -+ { -+ l_fact = L (); -+ a_fact = U (); // FIXME: sub-optimal -+ ipvt = getp (); -+ } -+} -+ -+template -+lu_type -+base_lu :: L (void) const -+{ -+ if (packed ()) -+ { -+ octave_idx_type a_nr = a_fact.rows (); -+ octave_idx_type a_nc = a_fact.cols (); -+ octave_idx_type mn = (a_nr < a_nc ? a_nr : a_nc); -+ -+ lu_type l (a_nr, mn, lu_elt_type (0.0)); -+ -+ for (octave_idx_type i = 0; i < a_nr; i++) -+ { -+ if (i < a_nc) -+ l.xelem (i, i) = 1.0; -+ -+ for (octave_idx_type j = 0; j < (i < a_nc ? i : a_nc); j++) -+ l.xelem (i, j) = a_fact.xelem (i, j); -+ } -+ -+ return l; -+ } -+ else -+ return l_fact; -+} -+ -+template -+lu_type -+base_lu :: U (void) const -+{ -+ if (packed ()) -+ { -+ octave_idx_type a_nr = a_fact.rows (); -+ octave_idx_type a_nc = a_fact.cols (); -+ octave_idx_type mn = (a_nr < a_nc ? a_nr : a_nc); -+ -+ lu_type u (mn, a_nc, lu_elt_type (0.0)); -+ -+ for (octave_idx_type i = 0; i < mn; i++) -+ { -+ for (octave_idx_type j = i; j < a_nc; j++) -+ u.xelem (i, j) = a_fact.xelem (i, j); -+ } -+ -+ return u; -+ } -+ else -+ return a_fact; -+} -+ -+template -+lu_type -+base_lu :: Y (void) const -+{ -+ if (! packed ()) -+ (*current_liboctave_error_handler) -+ ("lu: Y () not implemented for unpacked form"); -+ return a_fact; -+} -+ -+template -+Array -+base_lu :: getp (void) const -+{ -+ if (packed ()) -+ { -+ octave_idx_type a_nr = a_fact.rows (); -+ -+ Array pvt (dim_vector (a_nr, 1)); -+ -+ for (octave_idx_type i = 0; i < a_nr; i++) -+ pvt.xelem (i) = i; -+ -+ for (octave_idx_type i = 0; i < ipvt.length (); i++) -+ { -+ octave_idx_type k = ipvt.xelem (i); -+ -+ if (k != i) -+ { -+ octave_idx_type tmp = pvt.xelem (k); -+ pvt.xelem (k) = pvt.xelem (i); -+ pvt.xelem (i) = tmp; -+ } -+ } -+ -+ return pvt; -+ } -+ else -+ return ipvt; -+} -+ -+template -+PermMatrix -+base_lu :: P (void) const -+{ -+ return PermMatrix (getp (), false); -+} -+ -+template -+ColumnVector -+base_lu :: P_vec (void) const -+{ -+ octave_idx_type a_nr = a_fact.rows (); -+ -+ ColumnVector p (a_nr); -+ -+ Array pvt = getp (); -+ -+ for (octave_idx_type i = 0; i < a_nr; i++) -+ p.xelem (i) = static_cast (pvt.xelem (i) + 1); -+ -+ return p; -+} -+ -+template -+bool -+base_lu::regular (void) const -+{ -+ bool retval = true; -+ -+ octave_idx_type k = std::min (a_fact.rows (), a_fact.columns ()); -+ -+ for (octave_idx_type i = 0; i < k; i++) -+ { -+ if (a_fact(i, i) == lu_elt_type ()) -+ { -+ retval = false; -+ break; -+ } -+ } -+ -+ return retval; -+} -diff --git a/src/base-lu.h b/src/base-lu.h -new file mode 100644 ---- /dev/null -+++ b/src/base-lu.h -@@ -0,0 +1,87 @@ -+/* -+ -+Copyright (C) 1996-2015 John W. Eaton -+Copyright (C) 2009 VZLU Prague -+ -+This file is part of Octave. -+ -+Octave 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 3 of the License, or (at your -+option) any later version. -+ -+Octave 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 Octave; see the file COPYING. If not, see -+. -+ -+*/ -+ -+#if !defined (octave_base_lu_h) -+#define octave_base_lu_h 1 -+ -+#include "MArray.h" -+#include "dColVector.h" -+#include "PermMatrix.h" -+ -+template -+class -+base_lu -+{ -+public: -+ -+ typedef typename lu_type::element_type lu_elt_type; -+ -+ base_lu (void) -+ : a_fact (), l_fact (), ipvt () { } -+ -+ base_lu (const base_lu& a) -+ : a_fact (a.a_fact), l_fact (a.l_fact), ipvt (a.ipvt) { } -+ -+ base_lu (const lu_type& l, const lu_type& u, -+ const PermMatrix& p); -+ -+ base_lu& operator = (const base_lu& a) -+ { -+ if (this != &a) -+ { -+ a_fact = a.a_fact; -+ l_fact = a.l_fact; -+ ipvt = a.ipvt; -+ } -+ return *this; -+ } -+ -+ virtual ~base_lu (void) { } -+ -+ bool packed (void) const; -+ -+ void unpack (void); -+ -+ lu_type L (void) const; -+ -+ lu_type U (void) const; -+ -+ lu_type Y (void) const; -+ -+ PermMatrix P (void) const; -+ -+ ColumnVector P_vec (void) const; -+ -+ bool regular (void) const; -+ -+protected: -+ -+ Array getp (void) const; -+ -+ lu_type a_fact; -+ lu_type l_fact; -+ -+ Array ipvt; -+}; -+ -+#endif -diff --git a/src/galois-ops.h b/src/galois-ops.h ---- a/src/galois-ops.h -+++ b/src/galois-ops.h -@@ -21,6 +21,32 @@ - #if !defined (galois_octave_ops_h) - #define galois_octave_ops_h 1 - -+#if ! defined (CAST_BINOP_ARGS) -+# define CAST_BINOP_ARGS(t1, t2) \ -+ t1 v1 = dynamic_cast (a1); \ -+ t2 v2 = dynamic_cast (a2) -+#endif -+ -+#if ! defined (CAST_UNOP_ARG) -+# define CAST_UNOP_ARG(t) \ -+ t v = dynamic_cast (a) -+#endif -+ -+#if ! defined (BINOPDECL) -+# define BINOPDECL(name, a1, a2) \ -+ static octave_value \ -+ CONCAT2(oct_binop_, name) (const octave_base_value& a1, \ -+ const octave_base_value& a2) -+#endif -+ -+#if ! defined (CATOPDECL) -+# define CATOPDECL(name, a1, a2) \ -+ static octave_value \ -+ CONCAT2(oct_catop_, name) (octave_base_value& a1, \ -+ const octave_base_value& a2, \ -+ const Array& ra_idx) -+#endif -+ - // Override the operator and function definition defines from Octave - - #define DEFBINOP_OP_G(name, t1, t2, op) \ -diff --git a/src/galois.cc b/src/galois.cc ---- a/src/galois.cc -+++ b/src/galois.cc -@@ -27,7 +27,7 @@ - #include "galoisfield.h" - #include "galois-def.h" - --#include -+#include "base-lu.cc" - - galois_field_list stored_galois_fields; - -diff --git a/src/galois.h b/src/galois.h ---- a/src/galois.h -+++ b/src/galois.h -@@ -22,10 +22,10 @@ - #define octave_galois_int_h 1 - - #include --#include - #include - - #include "galoisfield.h" -+#include "base-lu.h" - - typedef void (*solve_singularity_handler) (double rcond); - -diff --git a/src/ov-galois.cc b/src/ov-galois.cc ---- a/src/ov-galois.cc -+++ b/src/ov-galois.cc -@@ -24,7 +24,6 @@ - #include - #include - #include --#include - #include - #include - #include -@@ -35,6 +34,22 @@ - #include "galois.h" - #include "ov-galois.h" - -+#if defined (HAVE_HDF5) -+# if defined (HAVE_HDF5_H) -+# include -+# endif -+# include "oct-hdf5-types.h" -+# if defined (OCTAVE_ENABLE_64) -+# define H5T_NATIVE_IDX H5T_NATIVE_INT64 -+# else -+# define H5T_NATIVE_IDX H5T_NATIVE_INT -+# endif -+#endif -+ -+#if ! defined (X_CAST) -+# define X_CAST(T, E) (T) (E) -+#endif -+ - #if defined (DEFINE_OCTAVE_ALLOCATOR) - DEFINE_OCTAVE_ALLOCATOR (octave_galois); - #endif -diff --git a/src/syndtable.cc b/src/syndtable.cc ---- a/src/syndtable.cc -+++ b/src/syndtable.cc -@@ -20,8 +20,8 @@ - - #include - --#define COL_MAJ(N) (N / (SIZEOF_INT << 3)) --#define COL_MIN(N) (N % (SIZEOF_INT << 3)) -+#define COL_MAJ(N) (N / (sizeof (int) << 3)) -+#define COL_MIN(N) (N % (sizeof (int) << 3)) - - Array - get_errs (const int& nmin, const int& nmax, const int &nerrs) diff -r 74fe6c154dee -r aad0cc165a27 src/of-communications-3-fixes.patch --- a/src/of-communications-3-fixes.patch Wed Jan 01 18:52:37 2020 -0500 +++ /dev/null Thu Jan 01 00:00:00 1970 +0000 @@ -1,153 +0,0 @@ -diff -uNr a/src/configure b/src/configure ---- a/src/configure 2015-04-04 12:29:43.708815162 -0400 -+++ b/src/configure 2018-01-09 16:18:21.287096634 -0500 -@@ -2450,7 +2450,6 @@ - cat confdefs.h - <<_ACEOF >conftest.$ac_ext - /* end confdefs.h. */ - -- #include - #include - - int -@@ -2502,7 +2501,6 @@ - cat confdefs.h - <<_ACEOF >conftest.$ac_ext - /* end confdefs.h. */ - -- #include - #include - class foo : public octave_base_value - { -@@ -2561,7 +2559,6 @@ - /* end confdefs.h. */ - - #include -- #include - #include - - int -diff -uNr a/src/configure.ac b/src/configure.ac ---- a/src/configure.ac 2018-01-09 16:19:20.488581117 -0500 -+++ b/src/configure.ac 2018-01-09 16:18:21.291096463 -0500 -@@ -34,7 +34,6 @@ - CXXFLAGS="$CXXFLAGS $comm_CXXFLAGS" - AC_LANG_PUSH(C++) - AC_COMPILE_IFELSE([AC_LANG_PROGRAM([[ -- #include - #include - ]], [[ - octave_hdf5_id x; -@@ -56,7 +55,6 @@ - CXXFLAGS="$CXXFLAGS $comm_CXXFLAGS" - AC_LANG_PUSH(C++) - AC_COMPILE_IFELSE([AC_LANG_PROGRAM([[ -- #include - #include - class foo : public octave_base_value - { -@@ -85,7 +83,6 @@ - AC_LANG_PUSH(C++) - AC_COMPILE_IFELSE([AC_LANG_PROGRAM([[ - #include -- #include - #include - ]], [[ - const octave_base_value x; x.print (std::cout); -diff -uNr a/src/galois.cc b/src/galois.cc ---- a/src/galois.cc 2018-01-09 16:19:22.260505820 -0500 -+++ b/src/galois.cc 2018-01-09 16:18:21.291096463 -0500 -@@ -18,7 +18,6 @@ - // program with any Open Source program, as defined by the Open Source - // Initiative (www.opensource.org) - --#include - #include - #include - #include -diff -uNr a/src/galois.h b/src/galois.h ---- a/src/galois.h 2018-01-09 16:19:22.260505820 -0500 -+++ b/src/galois.h 2018-01-09 16:18:46.190038466 -0500 -@@ -21,7 +21,6 @@ - #if !defined (octave_galois_int_h) - #define octave_galois_int_h 1 - --#include - #include - - #include "galoisfield.h" -diff -uNr a/src/gf.cc b/src/gf.cc ---- a/src/gf.cc 2015-04-04 12:28:43.954509931 -0400 -+++ b/src/gf.cc 2018-01-09 16:18:21.291096463 -0500 -@@ -29,7 +29,6 @@ - for a GPL release of his code - */ - --#include - #include - #include - #include -diff -uNr a/src/op-gm-gm.cc b/src/op-gm-gm.cc ---- a/src/op-gm-gm.cc 2015-04-04 12:28:43.954509931 -0400 -+++ b/src/op-gm-gm.cc 2018-01-09 16:18:21.295096293 -0500 -@@ -18,7 +18,6 @@ - // program with any Open Source program, as defined by the Open Source - // Initiative (www.opensource.org) - --#include - #include - #include - -diff -uNr a/src/op-gm-m.cc b/src/op-gm-m.cc ---- a/src/op-gm-m.cc 2015-04-04 12:28:43.954509931 -0400 -+++ b/src/op-gm-m.cc 2018-01-09 16:18:21.295096293 -0500 -@@ -18,7 +18,6 @@ - // program with any Open Source program, as defined by the Open Source - // Initiative (www.opensource.org) - --#include - #include - #include - -diff -uNr a/src/op-gm-s.cc b/src/op-gm-s.cc ---- a/src/op-gm-s.cc 2015-04-04 12:28:43.958509840 -0400 -+++ b/src/op-gm-s.cc 2018-01-09 16:18:21.295096293 -0500 -@@ -18,7 +18,6 @@ - // program with any Open Source program, as defined by the Open Source - // Initiative (www.opensource.org) - --#include - #include - #include - -diff -uNr a/src/op-m-gm.cc b/src/op-m-gm.cc ---- a/src/op-m-gm.cc 2015-04-04 12:28:43.958509840 -0400 -+++ b/src/op-m-gm.cc 2018-01-09 16:18:21.295096293 -0500 -@@ -18,7 +18,6 @@ - // program with any Open Source program, as defined by the Open Source - // Initiative (www.opensource.org) - --#include - #include - #include - -diff -uNr a/src/op-s-gm.cc b/src/op-s-gm.cc ---- a/src/op-s-gm.cc 2015-04-04 12:28:43.958509840 -0400 -+++ b/src/op-s-gm.cc 2018-01-09 16:18:21.295096293 -0500 -@@ -18,7 +18,6 @@ - // program with any Open Source program, as defined by the Open Source - // Initiative (www.opensource.org) - --#include - #include - #include - -diff -uNr a/src/ov-galois.cc b/src/ov-galois.cc ---- a/src/ov-galois.cc 2018-01-09 16:19:22.260505820 -0500 -+++ b/src/ov-galois.cc 2018-01-09 16:18:21.295096293 -0500 -@@ -20,7 +20,6 @@ - - #include - --#include - #include - #include - #include diff -r 74fe6c154dee -r aad0cc165a27 src/of-communications-4-fixes.patch --- a/src/of-communications-4-fixes.patch Wed Jan 01 18:52:37 2020 -0500 +++ /dev/null Thu Jan 01 00:00:00 1970 +0000 @@ -1,55 +0,0 @@ -diff -uNr a/src/base-lu.cc b/src/base-lu.cc ---- a/src/base-lu.cc 2018-01-09 18:06:56.797829797 -0500 -+++ b/src/base-lu.cc 2018-01-09 18:11:42.805590710 -0500 -@@ -21,10 +21,6 @@ - - */ - --#ifdef HAVE_CONFIG_H --#include --#endif -- - #include "base-lu.h" - - template -diff -uNr a/src/configure.ac b/src/configure.ac ---- a/src/configure.ac 2018-01-09 18:06:59.461715775 -0500 -+++ b/src/configure.ac 2018-01-09 18:25:15.510841692 -0500 -@@ -96,21 +96,5 @@ - PKG_CPPFLAGS="$PKG_CPPFLAGS -DHAVE_OCTAVE_BASE_VALUE_PRINT_CONST=1" - fi - --comm_save_CXXFLAGS=$CXXFLAGS --CXXFLAGS="$CXXFLAGS $comm_CXXFLAGS" --AC_LANG_PUSH(C++) --AC_CHECK_HEADER([octave/ls-oct-ascii.h], [PKG_CPPFLAGS="$PKG_CPPFLAGS -DHAVE_OCTAVE_ASCII_H=1"], [], --[ --#include --#incude --]) --AC_CHECK_HEADER([octave/ls-oct-text.h], [PKG_CPPFLAGS="$PKG_CPPFLAGS -DHAVE_OCTAVE_TEXT_H=1"], [], --[ --#include --#include --]) --AC_LANG_POP(C++) --CXXFLAGS=$comm_save_CXXFLAGS -- - AC_CONFIG_FILES([Makefile]) - AC_OUTPUT -diff -uNr a/src/ov-galois.cc b/src/ov-galois.cc ---- a/src/ov-galois.cc 2018-01-09 18:06:59.461715775 -0500 -+++ b/src/ov-galois.cc 2018-01-09 18:24:45.776112372 -0500 -@@ -29,11 +29,7 @@ - #include - - --#if defined(HAVE_OCTAVE_TEXT_H) --# include --#else --# include --#endif -+#include - - #include "galois.h" - #include "ov-galois.h" diff -r 74fe6c154dee -r aad0cc165a27 src/of-communications-5-fixes.patch --- a/src/of-communications-5-fixes.patch Wed Jan 01 18:52:37 2020 -0500 +++ /dev/null Thu Jan 01 00:00:00 1970 +0000 @@ -1,6045 +0,0 @@ -diff -uNr a/src/base-lu.cc b/src/base-lu.cc ---- a/src/base-lu.cc 2018-04-09 13:25:42.884981069 -0400 -+++ b/src/base-lu.cc 2018-04-09 13:53:52.501958203 -0400 -@@ -125,7 +125,7 @@ - for (octave_idx_type i = 0; i < a_nr; i++) - pvt.xelem (i) = i; - -- for (octave_idx_type i = 0; i < ipvt.length (); i++) -+ for (octave_idx_type i = 0; i < ipvt.numel (); i++) - { - octave_idx_type k = ipvt.xelem (i); - -diff -uNr a/src/base-lu.cc~ b/src/base-lu.cc~ ---- a/src/base-lu.cc~ 1969-12-31 19:00:00.000000000 -0500 -+++ b/src/base-lu.cc~ 2018-04-09 13:26:27.478896947 -0400 -@@ -0,0 +1,187 @@ -+/* -+ -+Copyright (C) 1996-2015 John W. Eaton -+Copyright (C) 2009 VZLU Prague -+ -+This file is part of Octave. -+ -+Octave 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 3 of the License, or (at your -+option) any later version. -+ -+Octave 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 Octave; see the file COPYING. If not, see -+. -+ -+*/ -+ -+#include "base-lu.h" -+ -+template -+base_lu::base_lu (const lu_type& l, const lu_type& u, -+ const PermMatrix& p) -+ : a_fact (u), l_fact (l), ipvt (p.transpose ().col_perm_vec ()) -+{ -+ if (l.columns () != u.rows ()) -+ (*current_liboctave_error_handler) ("lu: dimension mismatch"); -+} -+ -+template -+bool -+base_lu :: packed (void) const -+{ -+ return l_fact.dims () == dim_vector (); -+} -+ -+template -+void -+base_lu :: unpack (void) -+{ -+ if (packed ()) -+ { -+ l_fact = L (); -+ a_fact = U (); // FIXME: sub-optimal -+ ipvt = getp (); -+ } -+} -+ -+template -+lu_type -+base_lu :: L (void) const -+{ -+ if (packed ()) -+ { -+ octave_idx_type a_nr = a_fact.rows (); -+ octave_idx_type a_nc = a_fact.cols (); -+ octave_idx_type mn = (a_nr < a_nc ? a_nr : a_nc); -+ -+ lu_type l (a_nr, mn, lu_elt_type (0.0)); -+ -+ for (octave_idx_type i = 0; i < a_nr; i++) -+ { -+ if (i < a_nc) -+ l.xelem (i, i) = 1.0; -+ -+ for (octave_idx_type j = 0; j < (i < a_nc ? i : a_nc); j++) -+ l.xelem (i, j) = a_fact.xelem (i, j); -+ } -+ -+ return l; -+ } -+ else -+ return l_fact; -+} -+ -+template -+lu_type -+base_lu :: U (void) const -+{ -+ if (packed ()) -+ { -+ octave_idx_type a_nr = a_fact.rows (); -+ octave_idx_type a_nc = a_fact.cols (); -+ octave_idx_type mn = (a_nr < a_nc ? a_nr : a_nc); -+ -+ lu_type u (mn, a_nc, lu_elt_type (0.0)); -+ -+ for (octave_idx_type i = 0; i < mn; i++) -+ { -+ for (octave_idx_type j = i; j < a_nc; j++) -+ u.xelem (i, j) = a_fact.xelem (i, j); -+ } -+ -+ return u; -+ } -+ else -+ return a_fact; -+} -+ -+template -+lu_type -+base_lu :: Y (void) const -+{ -+ if (! packed ()) -+ (*current_liboctave_error_handler) -+ ("lu: Y () not implemented for unpacked form"); -+ return a_fact; -+} -+ -+template -+Array -+base_lu :: getp (void) const -+{ -+ if (packed ()) -+ { -+ octave_idx_type a_nr = a_fact.rows (); -+ -+ Array pvt (dim_vector (a_nr, 1)); -+ -+ for (octave_idx_type i = 0; i < a_nr; i++) -+ pvt.xelem (i) = i; -+ -+ for (octave_idx_type i = 0; i < ipvt.length (); i++) -+ { -+ octave_idx_type k = ipvt.xelem (i); -+ -+ if (k != i) -+ { -+ octave_idx_type tmp = pvt.xelem (k); -+ pvt.xelem (k) = pvt.xelem (i); -+ pvt.xelem (i) = tmp; -+ } -+ } -+ -+ return pvt; -+ } -+ else -+ return ipvt; -+} -+ -+template -+PermMatrix -+base_lu :: P (void) const -+{ -+ return PermMatrix (getp (), false); -+} -+ -+template -+ColumnVector -+base_lu :: P_vec (void) const -+{ -+ octave_idx_type a_nr = a_fact.rows (); -+ -+ ColumnVector p (a_nr); -+ -+ Array pvt = getp (); -+ -+ for (octave_idx_type i = 0; i < a_nr; i++) -+ p.xelem (i) = static_cast (pvt.xelem (i) + 1); -+ -+ return p; -+} -+ -+template -+bool -+base_lu::regular (void) const -+{ -+ bool retval = true; -+ -+ octave_idx_type k = std::min (a_fact.rows (), a_fact.columns ()); -+ -+ for (octave_idx_type i = 0; i < k; i++) -+ { -+ if (a_fact(i, i) == lu_elt_type ()) -+ { -+ retval = false; -+ break; -+ } -+ } -+ -+ return retval; -+} -diff -uNr a/src/cyclgen.cc b/src/cyclgen.cc ---- a/src/cyclgen.cc 2015-04-04 12:28:43.942510204 -0400 -+++ b/src/cyclgen.cc 2018-04-09 13:51:41.852070262 -0400 -@@ -29,7 +29,7 @@ - const Array& x, const int& n) - { - -- int x_len = x.length (); -+ int x_len = x.numel (); - Array si (dim_vector (n, 1), 0); - Array y (dim_vector (x_len, 1), 0); - -diff -uNr a/src/cyclgen.cc~ b/src/cyclgen.cc~ ---- a/src/cyclgen.cc~ 1969-12-31 19:00:00.000000000 -0500 -+++ b/src/cyclgen.cc~ 2015-04-04 12:28:43.942510204 -0400 -@@ -0,0 +1,278 @@ -+//Copyright (C) 2003 David Bateman -+// -+// 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 3 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; if not, see -+// . -+// -+// In addition to the terms of the GPL, you are permitted to link this -+// program with any Open Source program, as defined by the Open Source -+// Initiative (www.opensource.org) -+ -+#include -+ -+#include -+ -+// A simplified version of the filter function for specific lengths of a and b -+// in the Galois field GF(2) -+Array -+filter_gf2 (const Array& b, const Array& a, -+ const Array& x, const int& n) -+{ -+ -+ int x_len = x.length (); -+ Array si (dim_vector (n, 1), 0); -+ Array y (dim_vector (x_len, 1), 0); -+ -+ for (int i = 0; i < x_len; i++) -+ { -+ y(i) = si(0); -+ if (b(0) && x(i)) -+ y(i) ^= 1; -+ -+ for (int j = 0; j < n - 1; j++) -+ { -+ si(j) = si(j+1); -+ if (a(j+1) && y(i)) -+ si(j) ^= 1; -+ if (b(j+1) && x(i)) -+ si(j) ^= 1; -+ } -+ si(n-1) = 0; -+ if (a(n) && y(i)) -+ si(n-1) ^= 1; -+ if (b(n) && x(i)) -+ si(n-1) ^= 1; -+ } -+ -+ return y; -+} -+ -+// Cyclic polynomial is irreducible. I.E. it divides into x^n-1 -+// without remainder There must surely be an easier way of doing this -+// as the polynomials are over GF(2). -+static bool -+do_is_cyclic_polynomial (const Array& a, const int& n, const int& m) -+{ -+ Array y (dim_vector (n+1, 1), 0); -+ Array x (dim_vector (n-m+2, 1), 0); -+ y(0) = 1; -+ y(n) = 1; -+ x(0) = 1; -+ -+ Array b = filter_gf2 (y, a, x, n); -+ b.resize (dim_vector (n+1, 1), 0); -+ Array p (dim_vector (m+1, 1), 0); -+ p(0) = 1; -+ Array q = filter_gf2 (a, p, b, m); -+ -+ for (int i = 0; i < n+1; i++) -+ if (y(i) ^ q(i)) -+ return false; -+ -+ return true; -+} -+ -+DEFUN_DLD (cyclgen, args, nargout, -+ "-*- texinfo -*-\n\ -+@deftypefn {Loadable Function} {@var{h} =} cyclgen (@var{n}, @var{p})\n\ -+@deftypefnx {Loadable Function} {@var{h} =} cyclgen (@var{n}, @var{p}, @var{typ})\n\ -+@deftypefnx {Loadable Function} {[@var{h}, @var{g}] =} cyclgen (@dots{})\n\ -+@deftypefnx {Loadable Function} {[@var{h}, @var{g}, @var{k}] =} cyclgen (@dots{})\n\ -+Produce the parity check and generator matrix of a cyclic code. The parity\n\ -+check matrix is returned as a @var{m} by @var{n} matrix, representing the\n\ -+[@var{n},@var{k}] cyclic code. @var{m} is the order of the generator\n\ -+polynomial @var{p} and the message length @var{k} is given by\n\ -+@code{@var{n} - @var{m}}.\n\ -+\n\ -+The generator polynomial can either be a vector of ones and zeros,\n\ -+and length @var{m} representing,\n\ -+@tex\n\ -+$$ p_0 + p_1 x + p_2 x^2 + \\cdots + p_m x^{m-1} $$\n\ -+@end tex\n\ -+@ifnottex\n\ -+\n\ -+@example\n\ -+@var{p}(1) + @var{p}(2) * x + @var{p}(3) * x^2 + ... + @var{p}(@var{m}) * x^(m-1)\n\ -+@end example\n\ -+@end ifnottex\n\ -+\n\ -+The terms of the polynomial are stored least-significant term first.\n\ -+Alternatively, @var{p} can be an integer representation of the same\n\ -+polynomial.\n\ -+\n\ -+The form of the parity check matrix is determined by @var{typ}. If\n\ -+@var{typ} is 'system', a systematic parity check matrix is produced. If\n\ -+@var{typ} is 'nosys' and non-systematic parity check matrix is produced.\n\ -+\n\ -+If requested @code{cyclgen} also returns the @var{k} by @var{n} generator\n\ -+matrix @var{g}.\ -+\n\ -+@seealso{hammgen, gen2par, cyclpoly}\n\ -+@end deftypefn") -+{ -+ octave_value_list retval; -+ int nargin = args.length (); -+ unsigned long long p = 0; -+ int n, m, k, mm; -+ bool system = true; -+ Array pp; -+ -+ if (nargin < 2 || nargin > 3) -+ { -+ print_usage (); -+ return retval; -+ } -+ -+ n = args(0).int_value (); -+ m = 1; -+ while (n > (1<<(m+1))) -+ m++; -+ pp.resize (dim_vector (n+1, 1), 0); -+ -+ if (args(1).is_scalar_type ()) -+ { -+ p = (unsigned long long)(args(1).int_value ()); -+ mm = 1; -+ while (p > ((unsigned long long)1<<(mm+1))) -+ mm++; -+ for (int i = 0; i < mm+1; i++) -+ pp(i) = (p & (1< 2) -+ { -+ if (args(2).is_string ()) -+ { -+ std::string s_arg = args(2).string_value (); -+ -+ if (s_arg == "system") -+ system = true; -+ else if (s_arg == "nosys") -+ system = false; -+ else -+ { -+ error ("cyclgen: illegal argument"); -+ return retval; -+ } -+ } -+ else -+ { -+ error ("cyclgen: illegal argument"); -+ return retval; -+ } -+ } -+ -+ // Haven't implemented this since I'm not sure what matlab wants here -+ if (!system) -+ { -+ error ("cyclgen: non-systematic generator matrices not implemented"); -+ return retval; -+ } -+ -+ if (!do_is_cyclic_polynomial (pp, n, mm)) -+ { -+ error ("cyclgen: generator polynomial does not produce cyclic code"); -+ return retval; -+ } -+ -+ unsigned long long mask = 1; -+ unsigned long long *alpha_to = -+ (unsigned long long *)malloc (sizeof (unsigned long long) * n); -+ for (int i = 0; i < n; i++) -+ { -+ alpha_to[i] = mask; -+ mask <<= 1; -+ if (mask & ((unsigned long long)1< 1) -+ { -+ Matrix generator (k, n, 0); -+ -+ for (int i = 0; i < (int)k; i++) -+ for (int j = 0; j < (int)mm; j++) -+ generator(i, j) = parity(j, i+mm); -+ for (int i = 0; i < (int)k; i++) -+ generator(i, i+mm) = 1; -+ -+ retval(1) = octave_value (generator); -+ retval(2) = octave_value ((double)k); -+ } -+ return retval; -+} -+ -+/* -+%% Test input validation -+%!error cyclgen () -+%!error cyclgen (1) -+%!error cyclgen (1, 2, 3, 4) -+*/ -+ -+/* -+;;; Local Variables: *** -+;;; mode: C++ *** -+;;; End: *** -+*/ -diff -uNr a/src/cyclpoly.cc b/src/cyclpoly.cc ---- a/src/cyclpoly.cc 2015-04-04 12:28:43.942510204 -0400 -+++ b/src/cyclpoly.cc 2018-04-09 13:52:48.344959602 -0400 -@@ -38,7 +38,7 @@ - const Array& x, const int& n) - { - -- int x_len = x.length (); -+ int x_len = x.numel (); - Array si (dim_vector (n, 1), 0); - Array y (dim_vector (x_len, 1), 0); - -@@ -217,8 +217,8 @@ - for (unsigned long long i = (1UL<. -+// -+// In addition to the terms of the GPL, you are permitted to link this -+// program with any Open Source program, as defined by the Open Source -+// Initiative (www.opensource.org) -+ -+#include -+#include -+ -+#include -+ -+enum cyclic_poly_type -+{ -+ CYCLIC_POLY_MIN=0, -+ CYCLIC_POLY_MAX, -+ CYCLIC_POLY_ALL, -+ CYCLIC_POLY_L -+}; -+ -+// A simplified version of the filter function for specific lengths of -+// a and b in the Galois field GF(2) -+Array -+filter_gf2 (const Array& b, const Array& a, -+ const Array& x, const int& n) -+{ -+ -+ int x_len = x.numel (); -+ Array si (dim_vector (n, 1), 0); -+ Array y (dim_vector (x_len, 1), 0); -+ -+ for (int i = 0; i < x_len; i++) -+ { -+ y(i) = si(0); -+ if (b(0) && x(i)) -+ y(i) ^= 1; -+ -+ for (int j = 0; j < n - 1; j++) -+ { -+ si(j) = si(j+1); -+ if (a(j+1) && y(i)) -+ si(j) ^= 1; -+ if (b(j+1) && x(i)) -+ si(j) ^= 1; -+ } -+ si(n-1) = 0; -+ if (a(n) && y(i)) -+ si(n-1) ^= 1; -+ if (b(n) && x(i)) -+ si(n-1) ^= 1; -+ } -+ -+ return y; -+} -+ -+// Cyclic polynomial is irreducible. I.E. it divides into x^n-1 -+// without remainder There must surely be an easier way of doing this -+// as the polynomials are over GF(2). -+static bool -+do_is_cyclic_polynomial (const unsigned long long& a1, const int& n, -+ const int& m) -+{ -+ Array a (dim_vector (n+1, 1), 0); -+ Array y (dim_vector (n+1, 1), 0); -+ Array x (dim_vector (n-m+2, 1), 0); -+ y(0) = 1; -+ y(n) = 1; -+ x(0) = 1; -+ for (int i=0; i < m+1; i++) -+ a(i) = (a1 & (1UL << i) ? 1 : 0); -+ -+ Array b = filter_gf2 (y, a, x, n); -+ b.resize(dim_vector (n+1, 1), 0); -+ Array p (dim_vector (m+1, 1), 0); -+ p(0) = 1; -+ Array q = filter_gf2 (a, p, b, m); -+ -+ for (int i=0; i < n+1; i++) -+ if (y(i) ^ q(i)) -+ return false; -+ -+ return true; -+} -+ -+DEFUN_DLD (cyclpoly, args, nargout, -+ "-*- texinfo -*-\n\ -+@deftypefn {Loadable Function} {@var{y} =} cyclpoly (@var{n}, @var{k})\n\ -+@deftypefnx {Loadable Function} {@var{y} =} cyclpoly (@var{n}, @var{k}, @var{opt})\n\ -+@deftypefnx {Loadable Function} {@var{y} =} cyclpoly (@var{n}, @var{k}, @var{opt}, @var{rep})\n\ -+This function returns the cyclic generator polynomials of the code\n\ -+[@var{n},@var{k}]. By default the polynomial with the smallest weight\n\ -+is returned. However this behavior can be overridden with the @var{opt}\n\ -+flag. Valid values of @var{opt} are:\n\ -+\n\ -+@table @asis\n\ -+@item @code{\"all\"}\n\ -+Returns all of the polynomials of the code [@var{n},@var{k}]\n\ -+@item @code{\"min\"}\n\ -+Returns the polynomial of minimum weight of the code [@var{n},@var{k}]\n\ -+@item @code{\"max\"}\n\ -+Returns the polynomial of the maximum weight of the code [@var{n},@var{k}]\n\ -+@item @var{l}\n\ -+Returns the polynomials having exactly the weight @var{l}\n\ -+@end table\n\ -+\n\ -+The polynomials are returns as row-vectors in the variable @var{y}. Each\n\ -+row of @var{y} represents a polynomial with the least-significant term\n\ -+first. The polynomials can be returned with an integer representation\n\ -+if @var{rep} is @code{\"integer\"}. The default behavior is given if @var{rep}\n\ -+is @code{\"polynomial\"}.\n\ -+@seealso{gf, isprimitive}\n\ -+@end deftypefn") -+{ -+ octave_value retval; -+ int nargin = args.length (); -+ bool polyrep = true; -+ enum cyclic_poly_type type = CYCLIC_POLY_MIN; -+ RowVector cyclic_polys; -+ int l=0; -+ -+ if (nargin < 2 || nargin > 4) -+ { -+ print_usage (); -+ return retval; -+ } -+ -+ int n = args(0).int_value (); -+ int k = args(1).int_value ();; -+ -+ if (n < 1) -+ { -+ error ("cyclpoly: n must be 1 or greater"); -+ return retval; -+ } -+ -+ if (n <= k) -+ { -+ error ("cyclpoly: k must be less than n"); -+ return retval; -+ } -+ -+ for (int i = 2; i < nargin; i++) -+ { -+ if (args(i).is_scalar_type ()) -+ { -+ l = args(i).int_value (); -+ type = CYCLIC_POLY_L; -+ } -+ else if (args(i).is_string ()) -+ { -+ std::string s_arg = args(i).string_value (); -+ -+ if (s_arg == "integer") -+ polyrep = false; -+ else if (s_arg == "polynomial") -+ polyrep = true; -+ else if (s_arg == "min") -+ type = CYCLIC_POLY_MIN; -+ else if (s_arg == "max") -+ type = CYCLIC_POLY_MAX; -+ else if (s_arg == "all") -+ type = CYCLIC_POLY_ALL; -+ else -+ { -+ error ("cyclpoly: invalid argument"); -+ return retval; -+ } -+ } -+ else -+ { -+ error ("cyclpoly: incorrect argument type"); -+ return retval; -+ } -+ } -+ -+ int m = n - k; -+ -+ // Matlab code seems to think that 1+x+x^3 is of larger weight than -+ // 1+x^2+x^3. So for matlab compatiability the list of polynomials -+ // should be reversed by replacing "i+=2" with "i-=2" and visa-versa. -+ // Thats not going to happen!!! -+ -+ switch (type) -+ { -+ case CYCLIC_POLY_MIN: -+ cyclic_polys.resize (1); -+ for (unsigned long long i = (1UL< (1UL< --#include -+#include - #include - - #include "galois.h" -@@ -215,7 +215,7 @@ - - if (nr != a_nr || nc != a_nc) - { -- gripe_nonconformant ("operator +=", nr, nc, a_nr, a_nc); -+ octave::err_nonconformant ("operator +=", nr, nc, a_nr, a_nc); - return *this; - } - -@@ -251,7 +251,7 @@ - - if (nr != a_nr || nc != a_nc) - { -- gripe_nonconformant ("operator -=", nr, nc, a_nr, a_nc); -+ octave::err_nonconformant ("operator -=", nr, nc, a_nr, a_nc); - return *this; - } - -@@ -517,7 +517,7 @@ - { - if (a_nr != b_nr || a_nc != b_nc) - { -- gripe_nonconformant ("operator .^", a_nr, a_nc, a_nr, a_nc); -+ octave::err_nonconformant ("operator .^", a_nr, a_nc, a_nr, a_nc); - return galois (); - } - -@@ -548,7 +548,7 @@ - - if (a_nr != b_nr || a_nc != b_nc) - { -- gripe_nonconformant ("operator .^", a_nr, a_nc, b_nr, b_nc); -+ octave::err_nonconformant ("operator .^", a_nr, a_nc, b_nr, b_nc); - return galois (); - } - -@@ -757,7 +757,7 @@ - return product (a, b); - else if (a_nc != b_nr) - { -- gripe_nonconformant ("operator *", a_nr, a_nc, b_nr, b_nc); -+ octave::err_nonconformant ("operator *", a_nr, a_nc, b_nr, b_nc); - return galois (); - } - else -@@ -1302,7 +1302,7 @@ - - // Apply row interchanges to the right hand sides. - //for (int j = 0; j < IP.length (); j++) -- for (int j = IP.length ()-1; j >= 0; j--) -+ for (int j = IP.numel ()-1; j >= 0; j--) - { - int piv = IP(j); - for (int i = 0; i < b_nc; i++) -@@ -1334,7 +1334,7 @@ - Array IP (fact.ipvt); - - // Apply row interchanges to the right hand sides. -- for (int j = 0; j < IP.length (); j++) -+ for (int j = 0; j < IP.numel (); j++) - { - int piv = IP(j); - for (int i = 0; i < b_nc; i++) -@@ -1419,7 +1419,7 @@ - int a_nr = a.rows (); - int b_nr = b.rows (); - -- gripe_nonconformant ("operator /", a_nr, a_nc, b_nr, b_nc); -+ octave::err_nonconformant ("operator /", a_nr, a_nc, b_nr, b_nc); - return galois (); - } - -@@ -1463,7 +1463,7 @@ - int a_nc = a.cols (); - int b_nc = b.cols (); - -- gripe_nonconformant ("operator \\", a_nr, a_nc, b_nr, b_nc); -+ octave::err_nonconformant ("operator \\", a_nr, a_nc, b_nr, b_nc); - return galois (); - } - -diff -uNr a/src/galois.cc~ b/src/galois.cc~ ---- a/src/galois.cc~ 1969-12-31 19:00:00.000000000 -0500 -+++ b/src/galois.cc~ 2018-04-09 13:37:49.607001400 -0400 -@@ -0,0 +1,1494 @@ -+//Copyright (C) 2003 David Bateman -+// -+// 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 3 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; if not, see -+// . -+// -+// In addition to the terms of the GPL, you are permitted to link this -+// program with any Open Source program, as defined by the Open Source -+// Initiative (www.opensource.org) -+ -+#include -+#include -+#include -+ -+#include "galois.h" -+#include "galoisfield.h" -+#include "galois-def.h" -+ -+#include "base-lu.cc" -+ -+galois_field_list stored_galois_fields; -+ -+// galois class -+ -+galois::galois (const Array& a, const int& _m, -+ const int& _primpoly) : MArray (a.dims ()), field (NULL) -+{ -+ int _n = (1<<_m) - 1; -+ -+ // Check the validity of the data in the matrix -+ for (int i = 0; i < rows (); i++) -+ { -+ for (int j = 0; j < columns (); j++) -+ { -+ if ((a(i, j) < 0) || (a(i, j) > _n)) -+ { -+ gripe_range_galois (_m); -+ return; -+ } -+ xelem(i, j) = (int)a(i, j); -+ } -+ } -+ -+ field = stored_galois_fields.create_galois_field (_m, _primpoly); -+} -+ -+galois::galois (const MArray& a, const int& _m, -+ const int& _primpoly) : MArray (a.dims ()), field (NULL) -+{ -+ int _n = (1<<_m) - 1; -+ -+ // Check the validity of the data in the matrix -+ for (int i = 0; i < rows (); i++) -+ { -+ for (int j = 0; j < columns (); j++) -+ { -+ if ((a(i, j) < 0) || (a(i, j) > _n)) -+ { -+ gripe_range_galois (_m); -+ return; -+ } -+ xelem(i, j) = (int)a(i, j); -+ } -+ } -+ -+ field = stored_galois_fields.create_galois_field (_m, _primpoly); -+} -+ -+galois::galois (const Matrix& a, const int& _m, -+ const int& _primpoly) : MArray (a.dims ()), field (NULL) -+{ -+ int _n = (1<<_m) - 1; -+ -+ // Check the validity of the data in the matrix -+ for (int i = 0; i < rows (); i++) -+ { -+ for (int j = 0; j < columns (); j++) -+ { -+ if ((a(i, j) < 0) || (a(i, j) > _n)) -+ { -+ gripe_range_galois (_m); -+ return; -+ } -+ if ((a(i, j) - (double)((int)a(i, j))) != 0.) -+ { -+ gripe_integer_galois (); -+ return; -+ } -+ xelem(i, j) = (int)a(i, j); -+ } -+ } -+ -+ field = stored_galois_fields.create_galois_field (_m, _primpoly); -+} -+ -+galois::galois (int nr, int nc, const int& val, const int& _m, -+ const int& _primpoly) -+ : MArray (dim_vector (nr, nc), val), field (NULL) -+{ -+ int _n = (1<<_m) - 1; -+ -+ // Check the validity of the data in the matrix -+ if ((val < 0) || (val > _n)) -+ { -+ gripe_range_galois (_m); -+ return; -+ } -+ -+ field = stored_galois_fields.create_galois_field (_m, _primpoly); -+} -+ -+galois::galois (int nr, int nc, double val, const int& _m, -+ const int& _primpoly) -+ : MArray (dim_vector (nr, nc), (int)val), field (NULL) -+{ -+ int _n = (1<<_m) - 1; -+ -+ // Check the validity of the data in the matrix -+ if ((val < 0) || (val > _n)) -+ { -+ gripe_range_galois (_m); -+ return; -+ } -+ -+ if ((val - (double)((int)val)) != 0.) -+ { -+ gripe_integer_galois (); -+ return; -+ } -+ -+ field = stored_galois_fields.create_galois_field (_m, _primpoly); -+} -+ -+galois::galois (const galois& a) : MArray (a) -+{ -+ -+ if (!a.have_field ()) -+ { -+ gripe_copy_invalid_galois (); -+ field = NULL; -+ return; -+ } -+ -+ // This call to create_galois_field will just increment the usage counter -+ field = stored_galois_fields.create_galois_field (a.m (), a.primpoly ()); -+} -+ -+galois::~galois (void) -+{ -+ stored_galois_fields.delete_galois_field (field); -+ field = NULL; -+} -+ -+galois& -+galois::operator = (const galois& t) -+{ -+ if (!t.have_field ()) -+ { -+ gripe_copy_invalid_galois (); -+ if (have_field ()) -+ stored_galois_fields.delete_galois_field (field); -+ field = NULL; -+ return *this; -+ } -+ -+ if (have_field ()) -+ { -+ if ((m () != t.m ()) || (primpoly () != t.primpoly ())) -+ { -+ stored_galois_fields.delete_galois_field (field); -+ field = stored_galois_fields.create_galois_field (t.m (), t.primpoly ()); -+ } -+ } -+ else -+ field = stored_galois_fields.create_galois_field (t.m (), t.primpoly ()); -+ -+ // Copy the data -+ MArray::operator = (t); -+ -+ return *this; -+} -+ -+galois& -+galois::operator += (const galois& a) -+{ -+ int nr = rows (); -+ int nc = cols (); -+ -+ int a_nr = a.rows (); -+ int a_nc = a.cols (); -+ -+ if (have_field () && a.have_field ()) -+ { -+ if ((m () != a.m ()) || (primpoly () != a.primpoly ())) -+ { -+ gripe_differ_galois (); -+ return *this; -+ } -+ } -+ else -+ { -+ gripe_invalid_galois (); -+ return *this; -+ } -+ -+ if (nr != a_nr || nc != a_nc) -+ { -+ octave::err_nonconformant ("operator +=", nr, nc, a_nr, a_nc); -+ return *this; -+ } -+ -+ for (int i = 0; i < rows (); i++) -+ for (int j = 0; j < columns (); j++) -+ xelem(i, j) ^= a (i, j); -+ -+ return *this; -+} -+ -+galois& -+galois::operator -= (const galois& a) -+{ -+ int nr = rows (); -+ int nc = cols (); -+ -+ int a_nr = a.rows (); -+ int a_nc = a.cols (); -+ -+ if (have_field () && a.have_field ()) -+ { -+ if ((m () != a.m ()) || (primpoly () != a.primpoly ())) -+ { -+ gripe_differ_galois (); -+ return *this; -+ } -+ } -+ else -+ { -+ gripe_invalid_galois (); -+ return *this; -+ } -+ -+ if (nr != a_nr || nc != a_nc) -+ { -+ octave::err_nonconformant ("operator -=", nr, nc, a_nr, a_nc); -+ return *this; -+ } -+ -+ for (int i = 0; i < rows (); i++) -+ for (int j = 0; j < columns (); j++) -+ xelem(i, j) ^= a (i, j); -+ -+ return *this; -+} -+ -+galois -+galois::index (idx_vector& i, int resize_ok, const int& rfv) const -+{ -+ galois retval (MArray::index(i, resize_ok, rfv), m (), primpoly ()); -+ -+ return retval; -+} -+ -+galois -+galois::index (idx_vector& i, idx_vector& j, int resize_ok, -+ const int& rfv) const -+{ -+ galois retval (MArray::index(i, j, resize_ok, rfv), m (), primpoly ()); -+ -+ return retval; -+} -+ -+galois -+galois::concat (const galois& rb, const Array& ra_idx) -+{ -+ if (rb.numel () > 0) -+ insert (rb, ra_idx(0), ra_idx(1)); -+ return *this; -+} -+ -+galois -+galois::concat (const Matrix& rb, const Array& ra_idx) -+{ -+ if (numel () == 1) -+ return *this; -+ -+ galois tmp (0, 0, 0, m (), primpoly ()); -+ int _n = (1< _n)) -+ { -+ gripe_range_galois (m ()); -+ return *this; -+ } -+ if ((rb(i, j) - (double)((int)rb(i, j))) != 0.) -+ { -+ gripe_integer_galois (); -+ return *this; -+ } -+ tmp(i, j) = (int)rb(i, j); -+ } -+ } -+ -+ insert (tmp, ra_idx(0), ra_idx(1)); -+ return *this; -+} -+ -+galois -+concat (const Matrix& ra, const galois& rb, const Array& ra_idx) -+{ -+ galois retval (0, 0, 0, rb.m (), rb.primpoly ()); -+ int _n = (1< _n)) -+ { -+ gripe_range_galois (rb.m ()); -+ return retval; -+ } -+ if ((ra(i, j) - (double)((int)ra(i, j))) != 0.) -+ { -+ gripe_integer_galois (); -+ return retval; -+ } -+ retval(i, j) = (int)ra(i, j); -+#else -+ int tmp = (int)ra(i, j); -+ if (tmp < 0) -+ retval(i, j) = 0; -+ else if (tmp > _n) -+ retval(i, j) = _n; -+ else -+ retval(i, j) = tmp; -+#endif -+ } -+ } -+ -+ retval.insert (rb, ra_idx(0), ra_idx(1)); -+ return retval; -+} -+ -+galois& -+galois::insert (const galois& t, int r, int c) -+{ -+ if ((m () != t.m ()) || (primpoly () != t.primpoly ())) -+ (*current_liboctave_error_handler) ("inserted galois variable must " -+ "be in the same field"); -+ else -+ Array::insert (t, r, c); -+ return *this; -+} -+ -+galois -+galois::diag (void) const -+{ -+ return diag (0); -+} -+ -+galois -+galois::diag (int k) const -+{ -+ int nnr = rows (); -+ int nnc = cols (); -+ galois retval (0, 0, 0, m (), primpoly ()); -+ -+ if (k > 0) -+ nnc -= k; -+ else if (k < 0) -+ nnr += k; -+ -+ if (nnr > 0 && nnc > 0) -+ { -+ int ndiag = (nnr < nnc) ? nnr : nnc; -+ retval.resize (dim_vector (ndiag, 1)); -+ -+ if (k > 0) -+ { -+ for (int i = 0; i < ndiag; i++) -+ retval(i, 0) = xelem (i, i+k); -+ } -+ else if ( k < 0) -+ { -+ for (int i = 0; i < ndiag; i++) -+ retval(i, 0) = xelem (i-k, i); -+ } -+ else -+ { -+ for (int i = 0; i < ndiag; i++) -+ retval(i, 0) = xelem (i, i); -+ } -+ } -+ else -+ error ("diag: requested diagonal out of range"); -+ -+ return retval; -+} -+ -+// unary operations -+ -+boolMatrix -+galois::operator ! (void) const -+{ -+ int nr = rows (); -+ int nc = cols (); -+ -+ boolMatrix b (nr, nc); -+ -+ for (int j = 0; j < nc; j++) -+ for (int i = 0; i < nr; i++) -+ b (i, j) = ! xelem (i, j); -+ -+ return b; -+} -+ -+galois -+galois::transpose (void) const -+{ -+ galois a (Matrix (0, 0), m (), primpoly ()); -+ int d1 = rows (); -+ int d2 = cols (); -+ -+ a.resize (dim_vector (d2, d1)); -+ for (int j = 0; j < d2; j++) -+ for (int i = 0; i < d1; i++) -+ a (j, i) = xelem (i, j); -+ -+ return a; -+} -+ -+static inline int -+modn (int x, int m, int n) -+{ -+ while (x >= n) -+ { -+ x -= n; -+ x = (x >> m) + (x & n); -+ } -+ return x; -+} -+ -+galois -+elem_pow (const galois& a, const galois& b) -+{ -+ int a_nr = a.rows (); -+ int a_nc = a.cols (); -+ galois result (a_nr, a_nc, 0, a.m (), a.primpoly ()); -+ -+ int b_nr = b.rows (); -+ int b_nc = b.cols (); -+ -+ if (a.have_field () && b.have_field ()) -+ { -+ if ((a.m () != b.m ()) || (a.primpoly () != b.primpoly ())) -+ { -+ gripe_differ_galois (); -+ return galois (); -+ } -+ } -+ else -+ { -+ gripe_invalid_galois (); -+ return galois (); -+ } -+ -+ if (a_nr == 1 && a_nc == 1) -+ { -+ result.resize (dim_vector (b_nr, b_nc), 0); -+ int tmp = a.index_of (a(0, 0)); -+ for (int j = 0; j < b_nc; j++) -+ for (int i = 0; i < b_nr; i++) -+ if (b(i, j) == 0) -+ result(i, j) = 1; -+ else if (a(0, 0) != 0) -+ result(i, j) = a.alpha_to (modn (tmp * b(i, j), a.m (), a.n ())); -+ } -+ else if (b_nr == 1 && b_nc == 1) -+ { -+ for (int j = 0; j < a_nc; j++) -+ for (int i = 0; i < a_nr; i++) -+ if (b(0, 0) == 0) -+ result(i, j) = 1; -+ else if (a(i, j) != 0) -+ result(i, j) = a.alpha_to (modn (a.index_of (a(i, j)) * -+ b(0, 0), a.m (), a.n ())); -+ } -+ else -+ { -+ if (a_nr != b_nr || a_nc != b_nc) -+ { -+ octave::err_nonconformant ("operator .^", a_nr, a_nc, a_nr, a_nc); -+ return galois (); -+ } -+ -+ for (int j = 0; j < a_nc; j++) -+ for (int i = 0; i < a_nr; i++) -+ if (b(i, j) == 0) -+ result(i, j) = 1; -+ else if (a(i, j) != 0) -+ result(i, j) = a.alpha_to (modn (a.index_of (a(i, j)) * -+ b(i, j), a.m (), a.n ())); -+ } -+ -+ return result; -+} -+ -+galois -+elem_pow (const galois& a, const Matrix& b) -+{ -+ int a_nr = a.rows (); -+ int a_nc = a.cols (); -+ galois result (a_nr, a_nc, 0, a.m (), a.primpoly ()); -+ -+ int b_nr = b.rows (); -+ int b_nc = b.cols (); -+ -+ if (b_nr == 1 && b_nc == 1) -+ return elem_pow (a, b(0, 0)); -+ -+ if (a_nr != b_nr || a_nc != b_nc) -+ { -+ octave::err_nonconformant ("operator .^", a_nr, a_nc, b_nr, b_nc); -+ return galois (); -+ } -+ -+ for (int j = 0; j < a_nc; j++) -+ for (int i = 0; i < a_nr; i++) -+ { -+ int tmp = (int)b(i, j); -+ while (tmp < 0) -+ tmp += a.n (); -+ if (tmp == 0) -+ result(i, j) = 1; -+ else if (a(i, j) != 0) -+ result(i, j) = a.alpha_to (modn (a.index_of (a(i, j)) * tmp, -+ a.m (), a.n ())); -+ } -+ return result; -+} -+ -+galois -+elem_pow (const galois& a, double b) -+{ -+ int a_nr = a.rows (); -+ int a_nc = a.cols (); -+ galois result (a_nr, a_nc, 0, a.m (), a.primpoly ()); -+ int bi = (int) b; -+ -+ if ((double)bi != b) -+ { -+ gripe_integer_galois (); -+ return galois (); -+ } -+ -+ while (bi < 0) -+ bi += a.n (); -+ -+ for (int j = 0; j < a_nc; j++) -+ for (int i = 0; i < a_nr; i++) -+ { -+ if (bi == 0) -+ result(i, j) = 1; -+ else if (a(i, j) != 0) -+ result(i, j) = a.alpha_to (modn (a.index_of (a(i, j)) * -+ bi, a.m (), a.n ())); -+ } -+ return result; -+} -+ -+galois -+elem_pow (const galois &a, int b) -+{ -+ int a_nr = a.rows (); -+ int a_nc = a.cols (); -+ galois result (a_nr, a_nc, 0, a.m (), a.primpoly ()); -+ -+ while (b < 0) -+ b += a.n (); -+ -+ for (int j = 0; j < a_nc; j++) -+ for (int i = 0; i < a_nr; i++) -+ { -+ if (b == 0) -+ result(i, j) = 1; -+ else if (a(i, j) != 0) -+ result(i, j) = a.alpha_to (modn (a.index_of (a(i, j)) * b, -+ a.m (), a.n ())); -+ } -+ return result; -+} -+ -+galois -+pow (const galois& a, double b) -+{ -+ int bi = (int)b; -+ if ((double)bi != b) -+ { -+ gripe_integer_power_galois (); -+ return galois (); -+ } -+ -+ return pow (a, bi); -+} -+ -+galois -+pow (const galois& a, const galois& b) -+{ -+ int nr = b.rows (); -+ int nc = b.cols (); -+ -+ if (a.have_field () && b.have_field ()) -+ { -+ if ((a.m () != b.m ()) || (a.primpoly () != b.primpoly ())) -+ { -+ gripe_differ_galois (); -+ return galois (); -+ } -+ } -+ else -+ { -+ gripe_invalid_galois (); -+ return galois (); -+ } -+ -+ if (nr != 1 || nc != 1) -+ { -+ gripe_square_galois (); -+ return galois (); -+ } -+ else -+ return pow (a, b(0, 0)); -+} -+ -+galois -+pow (const galois& a, int b) -+{ -+ galois retval; -+ int nr = a.rows (); -+ int nc = a.cols (); -+ -+ if (!a.have_field ()) -+ { -+ gripe_invalid_galois (); -+ return retval; -+ } -+ -+ if (nr == 0 || nc == 0 || nr != nc) -+ gripe_square_galois (); -+ else if (b == 0) -+ { -+ retval = galois (nr, nc, 0, a.m (), a.primpoly ()); -+ for (int i = 0; i < nr; i++) -+ retval(i, i) = 1; -+ } -+ else -+ { -+ galois atmp; -+ -+ if (b < 0 ) -+ { -+ atmp = a.inverse (); -+ b = abs (b); -+ } -+ else -+ atmp = a; -+ -+ retval = atmp; -+ b--; -+ while (b > 0) -+ { -+ if (b & 1) -+ retval = retval * atmp; -+ -+ b >>= 1; -+ -+ if (b > 0) -+ atmp = atmp * atmp; -+ } -+ } -+ -+ return retval; -+} -+ -+galois -+operator * (const Matrix& a, const galois& b) -+{ -+ galois tmp (a, b.m (), b.primpoly ()); -+ -+ OCTAVE_QUIT; -+ -+ return tmp * b; -+} -+ -+galois -+operator * (const galois& a, const Matrix& b) -+{ -+ galois tmp (b, a.m (), a.primpoly ()); -+ -+ OCTAVE_QUIT; -+ -+ return a * tmp; -+} -+ -+galois -+operator * (const galois& a, const galois& b) -+{ -+ if (a.have_field () && b.have_field ()) -+ { -+ if ((a.m () != b.m ()) || (a.primpoly () != b.primpoly ())) -+ { -+ gripe_differ_galois (); -+ return galois (); -+ } -+ } -+ else -+ { -+ gripe_invalid_galois (); -+ return galois (); -+ } -+ -+ int a_nr = a.rows (); -+ int a_nc = a.cols (); -+ -+ int b_nr = b.rows (); -+ int b_nc = b.cols (); -+ -+ if ((a_nr == 1 && a_nc == 1) || (b_nr == 1 && b_nc == 1)) -+ return product (a, b); -+ else if (a_nc != b_nr) -+ { -+ octave::err_nonconformant ("operator *", a_nr, a_nc, b_nr, b_nc); -+ return galois (); -+ } -+ else -+ { -+ galois retval (a_nr, b_nc, 0, a.m (), a.primpoly ()); -+ if (a_nr != 0 && a_nc != 0 && b_nc != 0) -+ { -+ // This is not optimum for referencing b, but can use vector -+ // to represent index(a(k,j)). Seems to be the fastest. -+ galois c (a_nr, 1, 0, a.m (), a.primpoly ()); -+ for (int j = 0; j < b_nr; j++) -+ { -+ for (int k = 0; k < a_nr; k++) -+ c(k, 0) = a.index_of (a(k, j)); -+ -+ for (int i = 0; i < b_nc; i++) -+ if (b(j, i) != 0) -+ { -+ int tmp = a.index_of (b(j, i)); -+ for (int k = 0; k < a_nr; k++) -+ { -+ if (a(k, j) != 0) -+ retval(k, i) = retval(k, i) -+ ^ a.alpha_to (modn (tmp + c(k, 0), -+ a.m (), a.n ())); -+ } -+ } -+ } -+ } -+ return retval; -+ } -+} -+ -+// Other operators -+boolMatrix -+galois::all (int dim) const -+{ -+ return do_mx_red_op (*this, dim, mx_inline_all); -+} -+ -+boolMatrix -+galois::any (int dim) const -+{ -+ return do_mx_red_op (*this, dim, mx_inline_any); -+} -+ -+galois -+galois::prod (int dim) const -+{ -+ if (!have_field ()) -+ { -+ gripe_invalid_galois (); -+ return galois (); -+ } -+ -+ galois retval (0, 0, 0, m (), primpoly ()); -+ -+#define ROW_EXPR \ -+ if ((retval(i, 0) == 0) || (elem (i, j) == 0)) \ -+ retval(i, 0) = 0; \ -+ else \ -+ retval(i, 0) = alpha_to (modn (index_of (retval(i, 0)) + \ -+ index_of (elem (i, j)), m (), n ())); -+ -+#define COL_EXPR \ -+ if ((retval(0, j) == 0) || (elem (i, j) == 0)) \ -+ retval(0, j) = 0; \ -+ else \ -+ retval(0, j) = alpha_to (modn (index_of (retval(0, j)) + \ -+ index_of (elem (i, j)), m (), n ())); -+ -+ GALOIS_REDUCTION_OP (retval, ROW_EXPR, COL_EXPR, 1, 1); -+ return retval; -+ -+#undef ROW_EXPR -+#undef COL_EXPR -+} -+ -+galois -+galois::sum (int dim) const -+{ -+ if (!have_field ()) -+ { -+ gripe_invalid_galois (); -+ return galois (); -+ } -+ -+ galois retval (0, 0, 0, m (), primpoly ()); -+ -+ -+#define ROW_EXPR \ -+ retval(i, 0) ^= elem (i, j); -+ -+#define COL_EXPR \ -+ retval(0, j) ^= elem (i, j); -+ -+ GALOIS_REDUCTION_OP (retval, ROW_EXPR, COL_EXPR, 0, 0); -+ return retval; -+ -+#undef ROW_EXPR -+#undef COL_EXPR -+} -+ -+galois -+galois::sumsq (int dim) const -+{ -+ if (!have_field ()) -+ { -+ gripe_invalid_galois (); -+ return galois (); -+ } -+ -+ galois retval (0, 0, 0, m (), primpoly ()); -+ -+#define ROW_EXPR \ -+ if (elem (i, j) != 0) \ -+ retval(i, 0) ^= alpha_to (modn (2*index_of (elem (i, j)), m (), n ())); -+ -+#define COL_EXPR \ -+ if (elem (i, j) != 0) \ -+ retval(0, j) ^= alpha_to (modn (2*index_of (elem (i, j)), m (), n ())); -+ -+ GALOIS_REDUCTION_OP (retval, ROW_EXPR, COL_EXPR, 0, 0); -+ return retval; -+ -+#undef ROW_EXPR -+#undef COL_EXPR -+} -+ -+galois -+galois::sqrt (void) const -+{ -+ galois retval (*this); -+ int nr = rows (); -+ int nc = cols (); -+ -+ for (int j = 0; j < nc; j++) -+ { -+ for (int i = 0; i < nr; i++) -+ if (retval.index_of (retval(i, j)) & 1) -+ retval(i, j) = retval.alpha_to ((retval.index_of (retval(i, j)) -+ + retval.n ()) / 2); -+ else -+ retval(i, j) = retval.alpha_to (retval.index_of (retval(i, j)) -+ / 2); -+ } -+ return retval; -+} -+ -+galois -+galois::log (void) const -+{ -+ bool warned = false; -+ if (!have_field ()) -+ { -+ gripe_invalid_galois (); -+ return galois (); -+ } -+ -+ galois retval (*this); -+ int nr = rows (); -+ int nc = cols (); -+ -+ for (int j = 0; j < nc; j++) -+ for (int i = 0; i < nr; i++) -+ { -+ if (retval(i, j) == 0) -+ { -+ if (!warned) -+ { -+ warning ("log of zero undefined in Galois field"); -+ warned = true; -+ } -+ // How do I flag a NaN without either -+ // 1) Having to check everytime that the data is valid -+ // 2) Causing overflow in alpha_to or index_of!! -+ retval(i, j) = retval.index_of (retval(i, j)); -+ } -+ else -+ retval(i, j) = retval.index_of (retval(i, j)); -+ } -+ return retval; -+} -+ -+galois -+galois::exp (void) const -+{ -+ bool warned = false; -+ if (!have_field ()) -+ { -+ gripe_invalid_galois (); -+ return galois (); -+ } -+ -+ galois retval (*this); -+ int nr = rows (); -+ int nc = cols (); -+ -+ for (int j = 0; j < nc; j++) -+ for (int i = 0; i < nr; i++) -+ { -+ if (retval(i, j) == n ()) -+ { -+ if (!warned) -+ { -+ warning ("warning: exp of 2^m-1 undefined in Galois field"); -+ warned = true; -+ } -+ // How do I flag a NaN without either -+ // 1) Having to check everytime that the data is valid -+ // 2) Causing overflow in alpha_to or index_of!! -+ retval(i, j) = retval.alpha_to (retval(i, j)); -+ } -+ else -+ retval(i, j) = retval.alpha_to (retval(i, j)); -+ } -+ return retval; -+} -+ -+template class base_lu ; -+ -+void -+galoisLU::factor (const galois& a, const pivot_type& typ) -+{ -+ int a_nr = a.rows (); -+ int a_nc = a.cols (); -+ int mn = (a_nr > a_nc ? a_nc : a_nr); -+ -+ ptype = typ; -+ info = 0; -+ ipvt.resize (dim_vector (mn, 1)); -+ -+ a_fact = a; -+ -+ for (int j = 0; j < mn; j++) -+ { -+ int jp = j; -+ -+ // Find the pivot and test for singularity -+ if (ptype == galoisLU::ROW) -+ { -+ for (int i = j+1; i < a_nr; i++) -+ if (a_fact(i, j) > a_fact(jp, j)) -+ jp = i; -+ } -+ else -+ { -+ for (int i = j+1; i < a_nc; i++) -+ if (a_fact(j, i) > a_fact(j, jp)) -+ jp = i; -+ } -+ -+ ipvt(j) = jp; -+ -+ if (a_fact(jp, j) != 0) -+ { -+ if (ptype == galoisLU::ROW) -+ { -+ // Apply the interchange to columns 1:NC. -+ if (jp != j) -+ for (int i = 0; i < a_nc; i++) -+ { -+ int tmp = a_fact(j, i); -+ a_fact(j, i) = a_fact(jp, i); -+ a_fact(jp, i) = tmp; -+ } -+ } -+ else -+ { -+ // Apply the interchange to rows 1:NR. -+ if (jp != j) -+ for (int i = 0; i < a_nr; i++) -+ { -+ int tmp = a_fact(i, j); -+ a_fact(i, j) = a_fact(i, jp); -+ a_fact(i, jp) = tmp; -+ } -+ } -+ -+ // Compute elements J+1:M of J-th column. -+ if ( j < a_nr-1) -+ { -+ int idxj = a_fact.index_of (a_fact(j, j)); -+ for (int i = j+1; i < a_nr; i++) -+ { -+ if (a_fact(i, j) == 0) -+ a_fact(i, j) = 0; -+ else -+ a_fact(i, j) = a_fact.alpha_to (modn (a_fact.index_of (a_fact(i, j)) -+ - idxj + a_fact.n (), a_fact.m (), -+ a_fact.n ())); -+ } -+ } -+ } -+ else -+ { -+ info = 1; -+ } -+ -+ if (j < mn-1) -+ { -+ // Update trailing submatrix. -+ for (int i = j+1; i < a_nr; i++) -+ { -+ if (a_fact(i, j) != 0) -+ { -+ int idxi = a_fact.index_of (a_fact(i, j)); -+ for (int k = j+1; k < a_nc; k++) -+ { -+ if (a_fact(j, k) != 0) -+ a_fact(i, k) ^= a_fact.alpha_to (modn (a_fact.index_of (a_fact(j, k)) -+ + idxi, a_fact.m (), -+ a_fact.n ())); -+ } -+ } -+ } -+ } -+ } -+} -+ -+galois -+galoisLU::L (void) const -+{ -+ int a_nr = a_fact.rows (); -+ int a_nc = a_fact.cols (); -+ int mn = (a_nr < a_nc ? a_nr : a_nc); -+ -+ galois l (a_nr, mn, 0, a_fact.m (), a_fact.primpoly ()); -+ -+ for (int i = 0; i < mn; i++) -+ l(i, i) = 1; -+ -+ for (int j = 0; j < mn; j++) -+ for (int i = j+1; i < a_nr; i++) -+ l(i, j) = a_fact (i, j); -+ -+ return l; -+} -+ -+galois -+galoisLU::U (void) const -+{ -+ int a_nr = a_fact.rows (); -+ int a_nc = a_fact.cols (); -+ int mn = (a_nr < a_nc ? a_nr : a_nc); -+ -+ galois u (mn, a_nc, 0, a_fact.m (), a_fact.primpoly ()); -+ -+ for (int j = 0; j < a_nc; j++) -+ for (int i = 0; i < (j+1 > mn ? mn : j+1); i++) -+ u (i, j) = a_fact (i, j); -+ return u; -+} -+ -+galois -+galois::inverse (void) const -+{ -+ int info; -+ return inverse (info); -+} -+ -+galois -+galois::inverse (int& info, int force) const -+{ -+ int nr = rows (); -+ int nc = cols (); -+ info = 0; -+ -+ if (nr != nc || nr == 0 || nc == 0) -+ { -+ (*current_liboctave_error_handler) ("inverse requires square matrix"); -+ return galois (); -+ } -+ else -+ { -+ int info = 0; -+ -+ // Solve with identity matrix to find the inverse. -+ galois btmp (nr, nr, 0, m (), primpoly ()); -+ for (int i = 0; i < nr; i++) -+ btmp(i, i) = 1; -+ -+ galois retval = solve (btmp, info, 0); -+ -+ if (info == 0) -+ return retval; -+ else -+ return galois (); -+ } -+} -+ -+galois -+galois::determinant (void) const -+{ -+ int info; -+ return determinant (info); -+} -+ -+galois -+galois::determinant (int& info) const -+{ -+ galois retval (1, 1, 0, m (), primpoly ()); -+ -+ int nr = rows (); -+ int nc = cols (); -+ info = -1; -+ -+ if (nr == 0 || nc == 0) -+ { -+ info = 0; -+ retval(0, 0) = 1; -+ } -+ else -+ { -+ galoisLU fact (*this); -+ -+ if ( ! fact.singular ()) -+ { -+ galois A (fact.a_fact); -+ info = 0; -+ -+ retval(0, 0) = A(0, 0); -+ for (int i = 1; i < nr; i++) -+ { -+ if ((retval(0, 0) == 0) || (A(i, i) == 0)) -+ { -+ retval(0, 0) = 0; -+ error ("What the hell are we doing here!!!"); -+ } -+ else -+ retval(0, 0) = alpha_to (modn (index_of (retval(0, 0)) + -+ index_of (A(i, i)), m (), n ())); -+ } -+ } -+ } -+ -+ return retval; -+} -+ -+galois -+galois::solve (const galois& b) const -+{ -+ int info; -+ return solve (b, info); -+} -+ -+galois -+galois::solve (const galois& b, int& info) const -+{ -+ return solve (b, info, 0); -+} -+ -+galois -+galois::solve (const galois& b, int& info, -+ solve_singularity_handler sing_handler) const -+{ -+ galois retval (b); -+ -+ if (!have_field () || !b.have_field ()) -+ { -+ gripe_invalid_galois (); -+ return galois (); -+ } -+ else if ((m () != b.m ()) || (primpoly () != b.primpoly ())) -+ { -+ gripe_differ_galois (); -+ return galois (); -+ } -+ -+ int nr = rows (); -+ int nc = cols (); -+ int b_nr = b.rows (); -+ int b_nc = b.cols (); -+ galois c (nr, 1, 0, m (), primpoly ()); -+ -+ // if (nr == 0 || nc == 0 || nr != nc || nr != b_nr) -+ if (nr == 0 || nc == 0 || nr != b_nr) -+ { -+ (*current_liboctave_error_handler) -+ ("matrix dimension mismatch solution of linear equations"); -+ return galois (); -+ } -+ else if (nc > nr) -+ { -+ // Under-determined system, use column interchanges. -+ galoisLU fact ((*this), galoisLU::COL); -+ -+ if (fact.singular ()) -+ { -+ info = -1; -+ if (sing_handler) -+ sing_handler (0.0); -+ else -+ (*current_liboctave_error_handler)("galois matrix singular"); -+ -+ return galois (); -+ } -+ else -+ { -+ galois A (fact.a_fact); -+ Array IP (fact.ipvt); -+ -+ // Resize the number of solution rows if needed -+ if (nc > nr) -+ retval.resize (dim_vector (b_nr+nc-nr, b_nc), 0); -+ -+ //Solve L*X = B, overwriting B with X. -+ int mn = (nc < nr ? nc : nr); -+ for (int k = 0; k < mn; k++) -+ { -+ for (int i = k+1; i < nr; i++) -+ c(i, 0) = index_of (A(i, k)); -+ -+ for (int j = 0; j < b_nc; j++) -+ if (retval(k, j) != 0) -+ { -+ int idx = index_of (retval(k, j)); -+ for (int i = k+1; i < nr; i++) -+ if (A(i, k) != 0) -+ retval(i, j) ^= alpha_to (modn (c(i, 0) + idx, m (), n ())); -+ } -+ } -+ -+ // Solve U*X = B, overwriting B with X. -+ for (int k = (nc < nr ? nc-1 : nr-1); k >= 0; k--) -+ { -+ int mn = k+1 < nr ? k+1 : nr; -+ for (int i = 0; i < mn; i++) -+ c(i, 0) = index_of (A(i, k)); -+ mn = k < nr ? k : nr; -+ for (int j = 0; j < b_nc; j++) -+ if (retval(k, j) != 0) -+ { -+ retval(k, j) = alpha_to (modn (index_of (retval(k, j)) - -+ c(k, 0) + n (), m (), n ())); -+ int idx = index_of (retval(k, j)); -+ for (int i = 0; i < mn; i++) -+ if (A(i, k) != 0) -+ retval(i, j) ^= alpha_to (modn (c(i, 0) + idx, m (), n ())); -+ } -+ } -+ -+ // Apply row interchanges to the right hand sides. -+ //for (int j = 0; j < IP.length (); j++) -+ for (int j = IP.length ()-1; j >= 0; j--) -+ { -+ int piv = IP(j); -+ for (int i = 0; i < b_nc; i++) -+ { -+ int tmp = retval(j, i); -+ retval(j, i) = retval(piv, i); -+ retval(piv, i) = tmp; -+ } -+ } -+ } -+ } -+ else -+ { -+ galoisLU fact (*this); -+ -+ if (fact.singular ()) -+ { -+ info = -1; -+ if (sing_handler) -+ sing_handler (0.0); -+ else -+ (*current_liboctave_error_handler)("galois matrix singular"); -+ -+ return galois (); -+ } -+ else -+ { -+ galois A (fact.a_fact); -+ Array IP (fact.ipvt); -+ -+ // Apply row interchanges to the right hand sides. -+ for (int j = 0; j < IP.length (); j++) -+ { -+ int piv = IP(j); -+ for (int i = 0; i < b_nc; i++) -+ { -+ int tmp = retval(j, i); -+ retval(j, i) = retval(piv, i); -+ retval(piv, i) = tmp; -+ } -+ } -+ -+ //Solve L*X = B, overwriting B with X. -+ int mn = (nc < nr ? nc : nr); -+ for (int k = 0; k < mn; k++) -+ { -+ for (int i = k+1; i < nr; i++) -+ c(i, 0) = index_of (A(i, k)); -+ for (int j = 0; j < b_nc; j++) -+ if (retval(k, j) != 0) -+ { -+ int idx = index_of (retval(k, j)); -+ for (int i = k+1; i < nr; i++) -+ if (A(i, k) != 0) -+ retval(i, j) ^= alpha_to (modn (c(i, 0) + idx, m (), n ())); -+ } -+ } -+ -+ // Solve U*X = B, overwriting B with X. -+ for (int k = (nc < nr ? nc-1 : nr-1); k >= 0; k--) -+ { -+ int mn = k+1 < nr ? k+1 : nr; -+ for (int i = 0; i < mn; i++) -+ c(i, 0) = index_of (A(i, k)); -+ mn = k < nr ? k : nr; -+ for (int j = 0; j < b_nc; j++) -+ if (retval(k, j) != 0) -+ { -+ retval(k, j) = alpha_to (modn (index_of (retval(k, j)) - -+ c(k, 0) + n (), m (), n ())); -+ int idx = index_of (retval(k, j)); -+ for (int i = 0; i < mn; i++) -+ if (A(i, k) != 0) -+ retval(i, j) ^= alpha_to (modn (c(i, 0) + idx, m (), n ())); -+ } -+ } -+ -+ // Resize the number of solution rows if needed -+ if (nc < nr) -+ retval.resize (dim_vector (b_nr+nc-nr, b_nc)); -+ -+ } -+ } -+ -+ return retval; -+} -+ -+galois -+xdiv (const galois& a, const Matrix& b) -+{ -+ galois btmp (b, a.m (), a.primpoly ()); -+ -+ return xdiv (a, btmp); -+} -+ -+galois -+xdiv (const Matrix& a, const galois& b) -+{ -+ galois atmp (a, b.m (), b.primpoly ()); -+ -+ return xdiv (atmp, b); -+} -+ -+galois -+xdiv (const galois& a, const galois& b) -+{ -+ int info = 0; -+ int a_nc = a.cols (); -+ int b_nc = b.cols (); -+ -+ // if ((a_nc != b_nc) || (b.rows () != b.cols ())) -+ if (a_nc != b_nc) -+ { -+ int a_nr = a.rows (); -+ int b_nr = b.rows (); -+ -+ octave::err_nonconformant ("operator /", a_nr, a_nc, b_nr, b_nc); -+ return galois (); -+ } -+ -+ galois atmp = a.transpose (); -+ galois btmp = b.transpose (); -+ galois result = btmp.solve (atmp, info, 0); -+ -+ if (info == 0) -+ return galois (result.transpose ()); -+ else -+ return galois (); -+} -+ -+ -+galois -+xleftdiv (const galois& a, const Matrix& b) -+{ -+ galois btmp (b, a.m (), a.primpoly ()); -+ -+ return xleftdiv (a, btmp); -+} -+ -+galois -+xleftdiv (const Matrix& a, const galois& b) -+{ -+ galois atmp (a, b.m (), b.primpoly ()); -+ -+ return xleftdiv (atmp, b); -+} -+ -+galois -+xleftdiv (const galois& a, const galois& b) -+{ -+ int info = 0; -+ int a_nr = a.rows (); -+ int b_nr = b.rows (); -+ -+ // if ((a_nr != b_nr) || (a.rows () != a.columns ())) -+ if (a_nr != b_nr) -+ { -+ int a_nc = a.cols (); -+ int b_nc = b.cols (); -+ -+ octave::err_nonconformant ("operator \\", a_nr, a_nc, b_nr, b_nc); -+ return galois (); -+ } -+ -+ galois result = a.solve (b, info, 0); -+ -+ if (info == 0) -+ return result; -+ else -+ return galois (); -+} -+ -+MM_BIN_OPS1 (galois, galois, galois, 1, 2, GALOIS) -+MM_BIN_OPS1 (galois, galois, Matrix, 1, 2, MATRIX) -+MM_BIN_OPS1 (galois, Matrix, galois, 2, 1, MATRIX) -+ -+MM_CMP_OPS1 (galois, , galois, , 1, 2, GALOIS) -+MM_CMP_OPS1 (galois, , Matrix, , 1, 2, MATRIX) -+MM_CMP_OPS1 (Matrix, , galois, , 2, 1, MATRIX) -+ -+MM_BOOL_OPS1 (galois, galois, 0.0, 1, 2, GALOIS) -+MM_BOOL_OPS1 (galois, Matrix, 0.0, 1, 2, MATRIX) -+MM_BOOL_OPS1 (Matrix, galois, 0.0, 2, 1, MATRIX) -+ -+/* -+;;; Local Variables: *** -+;;; mode: C++ *** -+;;; End: *** -+*/ -diff -uNr a/src/galois-def.h b/src/galois-def.h ---- a/src/galois-def.h 2015-04-04 12:28:43.942510204 -0400 -+++ b/src/galois-def.h 2018-04-09 13:37:40.547425139 -0400 -@@ -137,7 +137,7 @@ - r(i, j) = (int)m1(i, j) ^ (int)m2(0, 0); \ - } \ - else \ -- gripe_nonconformant (#OP, m1_nr, m1_nc, m2_nr, m2_nc); \ -+ octave::err_nonconformant (#OP, m1_nr, m1_nc, m2_nr, m2_nc); \ - } \ - else \ - { \ -@@ -221,7 +221,7 @@ - } \ - } \ - else \ -- gripe_nonconformant (#F, m1_nr, m1_nc, m2_nr, m2_nc); \ -+ octave::err_nonconformant (#F, m1_nr, m1_nc, m2_nr, m2_nc); \ - } \ - else \ - if (m1_nr > 0 && m1_nc > 0) \ -@@ -289,7 +289,7 @@ - r(i, j) = C1 (m1(i, j)) OP C2 (m2(0, 0)); \ - } \ - else \ -- gripe_nonconformant (#F, m1_nr, m1_nc, m2_nr, m2_nc); \ -+ octave::err_nonconformant (#F, m1_nr, m1_nc, m2_nr, m2_nc); \ - } \ - \ - return r; \ -@@ -350,7 +350,7 @@ - OP (m2(0, 0) != ZERO); \ - } \ - else if ((m1_nr != 0 || m1_nc != 0) && (m2_nr != 0 || m2_nc != 0)) \ -- gripe_nonconformant (#F, m1_nr, m1_nc, m2_nr, m2_nc); \ -+ octave::err_nonconformant (#F, m1_nr, m1_nc, m2_nr, m2_nc); \ - } \ - \ - return r; \ -diff -uNr a/src/genqamdemod.cc b/src/genqamdemod.cc ---- a/src/genqamdemod.cc 2015-04-04 12:28:43.950510022 -0400 -+++ b/src/genqamdemod.cc 2018-04-09 13:33:02.852412423 -0400 -@@ -36,7 +36,7 @@ - - int nr1 (args(0).rows ()); - int nc1 (args(0).columns ()); -- int arg_is_empty1 = empty_arg ("genqamdemod", nr1, nc1); -+ int arg_is_empty1 = args(0).isempty (); - Matrix y (nr1,nc1); - - int nr2 (args(1).rows ()); -@@ -48,7 +48,7 @@ - if (arg_is_empty1 > 0) - return octave_value (Matrix ()); - -- if (args(0).is_real_type () && args(1).is_real_type ()) -+ if (args(0).isreal () && args(1).isreal ()) - { // Real-valued signal & constellation - Matrix x (args(0).matrix_value ()); - ColumnVector constellation (args(1).vector_value ()); -@@ -70,7 +70,7 @@ - } - } - } -- else if (args(0).is_complex_type () || args(1).is_complex_type ()) -+ else if (args(0).iscomplex () || args(1).iscomplex ()) - { // Complex-valued input & constellation - ComplexMatrix x (args(0).complex_matrix_value ()); - ComplexColumnVector constellation (args(1).complex_vector_value ()); -diff -uNr a/src/gf.cc b/src/gf.cc ---- a/src/gf.cc 2018-04-09 13:25:42.880981256 -0400 -+++ b/src/gf.cc 2018-04-09 14:16:54.029455163 -0400 -@@ -30,7 +30,8 @@ - */ - - #include --#include -+#include -+#include - #include - #include - #include -@@ -72,7 +73,7 @@ - // functions, this can't be done at the point. So if more default primitive - // polynomials are added to galoisfield.cc, need to update the "16" here - // as well!! --DEFUN_DLD (gf, args, nargout, -+DEFMETHOD_DLD (gf, interp, args, nargout, - "-*- texinfo -*-\n\ - @deftypefn {Loadable Function} {@var{y} =} gf (@var{x})\n\ - @deftypefnx {Loadable Function} {@var{y} =} gf (@var{x}, @var{m})\n\ -@@ -121,7 +122,7 @@ - install_s_gm_ops (); - install_gm_s_ops (); - galois_type_loaded = true; -- mlock (); -+ interp.mlock (); - } - - retval = new octave_galois (data, m, primpoly); -@@ -141,7 +142,7 @@ - - if ((!galois_type_loaded) || (a.type_id () != - octave_galois::static_type_id ())) -- gripe_wrong_type_arg ("gdiag", a); -+ err_wrong_type_arg ("gdiag", a); - else - { - galois m = ((const octave_galois&) a.get_rep ()).galois_value (); -@@ -190,12 +191,12 @@ - else - { - galois r = m.diag (k); -- if (r.capacity () > 0) -+ if (r.numel () > 0) - retval = new octave_galois (r); - } - } - else -- gripe_wrong_type_arg ("gdiag", a); -+ err_wrong_type_arg ("gdiag", a); - } - return retval; - } -@@ -301,7 +302,7 @@ - if ((!galois_type_loaded) || (args(0).type_id () != - octave_galois::static_type_id ())) - { -- gripe_wrong_type_arg ("greshape", args(0)); -+ err_wrong_type_arg ("greshape", args(0)); - return retval; - } - galois a = ((const octave_galois&) args(0).get_rep ()).galois_value (); -@@ -371,7 +372,7 @@ - } \ - else \ - { \ -- gripe_wrong_type_arg (#FCN, arg); \ -+ err_wrong_type_arg (#FCN, arg); \ - return retval; \ - } \ - } \ -@@ -468,7 +469,7 @@ - if (!galois_type_loaded || (args(0).type_id () != - octave_galois::static_type_id ())) - { -- gripe_wrong_type_arg ("gsqrt", args(0)); -+ err_wrong_type_arg ("gsqrt", args(0)); - return retval; - } - -@@ -507,7 +508,7 @@ - if (!galois_type_loaded || (args(0).type_id () != - octave_galois::static_type_id ())) - { -- gripe_wrong_type_arg ("glog", args(0)); -+ err_wrong_type_arg ("glog", args(0)); - return retval; - } - -@@ -546,7 +547,7 @@ - if (!galois_type_loaded || (args(0).type_id () != - octave_galois::static_type_id ())) - { -- gripe_wrong_type_arg ("gexp", args(0)); -+ err_wrong_type_arg ("gexp", args(0)); - return retval; - } - -@@ -577,9 +578,9 @@ - galois - filter (galois& b, galois& a, galois& x, galois& si) - { -- int ab_len = (a.length () > b.length () ? a.length () : b.length ()); -+ int ab_len = (a.numel () > b.numel () ? a.numel () : b.numel ()); - b.resize (dim_vector (ab_len, 1), 0); -- galois retval (x.length (), 1, 0, b.m (), b.primpoly ()); -+ galois retval (x.numel (), 1, 0, b.m (), b.primpoly ()); - int norm = a(0, 0); - - if (norm == 0) -@@ -587,43 +588,43 @@ - error ("gfilter: the first element of a must be non-zero"); - return galois (); - } -- if (si.length () != ab_len - 1) -+ if (si.numel () != ab_len - 1) - { -- error ("gfilter: si must be a vector of length max(length(a), length(b)) - 1"); -+ error ("gfilter: si must be a vector of length max(numel(a), numel(b)) - 1"); - return galois (); - } - if (norm != 1) - { - int idx_norm = b.index_of (norm); -- for (int i = 0; i < b.length (); i++) -+ for (int i = 0; i < b.numel (); i++) - { - if (b(i, 0) != 0) - b(i, 0) = b.alpha_to (modn (b.index_of (b(i, 0))-idx_norm+b.n (), - b.m (), b.n ())); - } - } -- if (a.length () > 1) -+ if (a.numel () > 1) - { - a.resize (dim_vector (ab_len, 1), 0); - - if (norm != 1) - { - int idx_norm = a.index_of (norm); -- for (int i = 0; i < a.length (); i++) -+ for (int i = 0; i < a.numel (); i++) - if (a(i, 0) != 0) - a(i, 0) = a.alpha_to (modn (a.index_of (a(i, 0))-idx_norm+a.n (), - a.m (), a.n ())); - } - -- for (int i = 0; i < x.length (); i++) -+ for (int i = 0; i < x.numel (); i++) - { - retval(i, 0) = si(0, 0); - if ((b(0, 0) != 0) && (x(i, 0) != 0)) - retval(i, 0) ^= b.alpha_to (modn (b.index_of (b(0, 0)) + - b.index_of (x(i, 0)), b.m (), b.n ())); -- if (si.length () > 1) -+ if (si.numel () > 1) - { -- for (int j = 0; j < si.length () - 1; j++) -+ for (int j = 0; j < si.numel () - 1; j++) - { - si(j, 0) = si(j+1, 0); - if ((a(j+1, 0) != 0) && (retval(i, 0) != 0)) -@@ -633,13 +634,13 @@ - si(j, 0) ^= b.alpha_to (modn (b.index_of (b(j+1, 0)) + - b.index_of (x(i, 0)), b.m (), b.n ())); - } -- si(si.length ()-1, 0) = 0; -- if ((a(si.length (), 0) != 0) && (retval(i, 0) != 0)) -- si(si.length ()-1, 0) ^= a.alpha_to (modn (a.index_of (a(si.length (), 0)) -+ si(si.numel ()-1, 0) = 0; -+ if ((a(si.numel (), 0) != 0) && (retval(i, 0) != 0)) -+ si(si.numel ()-1, 0) ^= a.alpha_to (modn (a.index_of (a(si.numel (), 0)) - + a.index_of (retval(i, 0)), - a.m (), a.n ())); -- if ((b(si.length (), 0) != 0) && (x(i, 0) != 0)) -- si(si.length ()-1, 0) ^= b.alpha_to (modn (b.index_of (b(si.length (), 0)) -+ if ((b(si.numel (), 0) != 0) && (x(i, 0) != 0)) -+ si(si.numel ()-1, 0) ^= b.alpha_to (modn (b.index_of (b(si.numel (), 0)) - + b.index_of (x(i, 0)), - b.m (), b.n ())); - } -@@ -655,26 +656,26 @@ - } - } - } -- else if (si.length () > 0) -+ else if (si.numel () > 0) - { -- for (int i = 0; i < x.length (); i++) -+ for (int i = 0; i < x.numel (); i++) - { - retval(i, 0) = si(0, 0); - if ((b(0, 0) != 0) && (x(i, 0) != 0)) - retval(i, 0) ^= b.alpha_to (modn (b.index_of (b(0, 0)) + - b.index_of (x(i, 0)), b.m (), b.n ())); -- if (si.length () > 1) -+ if (si.numel () > 1) - { -- for (int j = 0; j < si.length () - 1; j++) -+ for (int j = 0; j < si.numel () - 1; j++) - { - si(j, 0) = si(j+1, 0); - if ((b(j+1, 0) != 0) && (x(i, 0) != 0)) - si(j, 0) ^= b.alpha_to (modn (b.index_of (b(j+1, 0)) + - b.index_of (x(i, 0)), b.m (), b.n ())); - } -- si(si.length ()-1, 0) = 0; -- if ((b(si.length (), 0) != 0) && (x(i, 0) != 0)) -- si(si.length ()-1, 0) ^= b.alpha_to (modn (b.index_of (b(si.length (), 0)) -+ si(si.numel ()-1, 0) = 0; -+ if ((b(si.numel (), 0) != 0) && (x(i, 0) != 0)) -+ si(si.numel ()-1, 0) ^= b.alpha_to (modn (b.index_of (b(si.numel (), 0)) - + b.index_of (x(i, 0)), - b.m (), b.n ())); - } -@@ -688,7 +689,7 @@ - } - } - else -- for (int i = 0; i < x.length (); i++) -+ for (int i = 0; i < x.numel (); i++) - if ((b(0, 0) != 0) && (x(i, 0) != 0)) - retval(i, 0) = b.alpha_to (modn (b.index_of (b(0, 0)) + - b.index_of (x(i, 0)), b.m (), b.n ())); -@@ -717,7 +718,7 @@ - @smallexample\n\ - @group\n\ - N M\n\ -- SUM a(k+1) y(n-k) = SUM b(k+1) x(n-k) for 1<=n<=length(x)\n\ -+ SUM a(k+1) y(n-k) = SUM b(k+1) x(n-k) for 1<=n<=numel(x)\n\ - k=0 k=0\n\ - @end group\n\ - @end smallexample\n\ -@@ -729,7 +730,7 @@ - $a \\in \\Re^{N-1}$, $b \\in \\Re^{M-1}$, and $x \\in \\Re^P$.\n\ - @end tex\n\ - @ifnottex\n\ -- N=length(a)-1 and M=length(b)-1.\n\ -+ N=numel(a)-1 and M=numel(b)-1.\n\ - @end ifnottex\n\ - An equivalent form of this equation is:\n\ - @tex\n\ -@@ -743,7 +744,7 @@ - @smallexample\n\ - @group\n\ - N M\n\ -- y(n) = - SUM c(k+1) y(n-k) + SUM d(k+1) x(n-k) for 1<=n<=length(x)\n\ -+ y(n) = - SUM c(k+1) y(n-k) + SUM d(k+1) x(n-k) for 1<=n<=numel(x)\n\ - k=1 k=0\n\ - @end group\n\ - @end smallexample\n\ -@@ -838,8 +839,8 @@ - } - else - { -- int a_len = a.length (); -- int b_len = b.length (); -+ int a_len = a.numel (); -+ int b_len = b.numel (); - - int si_len = (a_len > b_len ? a_len : b_len) - 1; - -@@ -961,7 +962,7 @@ - if (!galois_type_loaded || (arg.type_id () != - octave_galois::static_type_id ())) - { -- gripe_wrong_type_arg ("glu", arg); -+ err_wrong_type_arg ("glu", arg); - return retval; - } - -@@ -970,7 +971,7 @@ - int nr = arg.rows (); - int nc = arg.columns (); - -- int arg_is_empty = empty_arg ("glu", nr, nc); -+ int arg_is_empty = arg.isempty (); - - if (arg_is_empty < 0) - return retval; -@@ -1048,13 +1049,13 @@ - if (!galois_type_loaded || (arg.type_id () != - octave_galois::static_type_id ())) - { -- gripe_wrong_type_arg ("ginverse", arg); -+ err_wrong_type_arg ("ginverse", arg); - return retval; - } - - galois m = ((const octave_galois&) arg.get_rep ()).galois_value (); - -- int arg_is_empty = empty_arg ("ginverse", nr, nc); -+ int arg_is_empty = arg.isempty (); - - if (arg_is_empty < 0) - return retval; -@@ -1065,7 +1066,7 @@ - } - if (nr != nc) - { -- gripe_square_matrix_required ("ginverse"); -+ err_square_matrix_required ("ginverse", "X"); - return retval; - } - -@@ -1142,7 +1143,7 @@ - if (!galois_type_loaded || (arg.type_id () != - octave_galois::static_type_id ())) - { -- gripe_wrong_type_arg ("gdet", arg); -+ err_wrong_type_arg ("gdet", arg); - return retval; - } - -@@ -1151,7 +1152,7 @@ - - galois m = ((const octave_galois&) arg.get_rep ()).galois_value (); - -- int arg_is_empty = empty_arg ("gdet", nr, nc); -+ int arg_is_empty = arg.isempty (); - - if (arg_is_empty < 0) - return retval; -@@ -1163,7 +1164,7 @@ - - if (nr != nc) - { -- gripe_square_matrix_required ("det"); -+ err_square_matrix_required ("det", "A"); - return retval; - } - -@@ -1202,7 +1203,7 @@ - if (!galois_type_loaded || (arg.type_id () != - octave_galois::static_type_id ())) - { -- gripe_wrong_type_arg ("grank", arg); -+ err_wrong_type_arg ("grank", arg); - return retval; - } - -@@ -1211,7 +1212,7 @@ - - galois m = ((const octave_galois&) arg.get_rep ()).galois_value (); - -- int arg_is_empty = empty_arg ("grank", nr, nc); -+ int arg_is_empty = arg.isempty (); - - if (arg_is_empty > 0) - retval = 0.0; -@@ -1332,7 +1333,7 @@ - if (!galois_type_loaded || (args(0).type_id () != - octave_galois::static_type_id ())) - { -- gripe_wrong_type_arg ("rsenc", args(0)); -+ err_wrong_type_arg ("rsenc", args(0)); - return retval; - } - -@@ -1873,7 +1874,7 @@ - if (!galois_type_loaded || (args(0).type_id () != - octave_galois::static_type_id ())) - { -- gripe_wrong_type_arg ("rsdec", args(0)); -+ err_wrong_type_arg ("rsdec", args(0)); - return retval; - } - -diff -uNr a/src/gf.cc~ b/src/gf.cc~ ---- a/src/gf.cc~ 1969-12-31 19:00:00.000000000 -0500 -+++ b/src/gf.cc~ 2018-04-09 14:12:50.168734516 -0400 -@@ -0,0 +1,2809 @@ -+// Copyright (C) 1994-1997 Robert Morelos-Zaragoza -+// Copyright (C) 2002 Phil Karn -+// Copyright (C) 2003 David Bateman -+// -+// 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 3 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; if not, see -+// . -+// -+// In addition to the terms of the GPL, you are permitted to link this -+// program with any Open Source program, as defined by the Open Source -+// Initiative (www.opensource.org) -+ -+/* -+Part of the function rsenc and the function decode_rs are from Phil Karn. See -+the website http://www.ka9q.net/code/fec for more details. -+ -+Parts of the function bchenco and bchdeco are from Robert Morelos-Zaragoza. See -+the website http://www.eccpage.com for more details. Permission has been granted -+for a GPL release of his code -+*/ -+ -+#include -+#include -+#include -+#include -+#include -+#include -+#include -+ -+#include "galois.h" -+#include "ov-galois.h" -+ -+static bool galois_type_loaded = false; -+ -+// PKG_ADD: autoload ("isgalois", "gf.oct"); -+// PKG_DEL: autoload ("isgalois", "gf.oct", "remove"); -+DEFUN_DLD (isgalois, args, , -+ "-*- texinfo -*-\n\ -+@deftypefn {Loadable Function} {} isgalois (@var{expr})\n\ -+Return 1 if the value of the expression @var{expr} is a Galois Field.\n\ -+@end deftypefn") -+{ -+ if (args.length () != 1) -+ print_usage (); -+ else if (!galois_type_loaded) -+ // Can be of Galois type if the type isn't load :-/ -+ return octave_value (0.); -+ else -+ return octave_value (args(0).type_id () == -+ octave_galois::static_type_id ()); -+ return octave_value (); -+} -+ -+/* -+%% Test input validation -+%!error isgalois () -+%!error isgalois (1, 2) -+*/ -+ -+// FIXME: -+// I want to replace the "16" below with __OCTAVE_GALOIS_MAX_M_AS_STRING, -+// but as I don't run the preprocessor when getting the help from the -+// functions, this can't be done at the point. So if more default primitive -+// polynomials are added to galoisfield.cc, need to update the "16" here -+// as well!! -+DEFMETHOD_DLD (interp, gf, args, nargout, -+ "-*- texinfo -*-\n\ -+@deftypefn {Loadable Function} {@var{y} =} gf (@var{x})\n\ -+@deftypefnx {Loadable Function} {@var{y} =} gf (@var{x}, @var{m})\n\ -+@deftypefnx {Loadable Function} {@var{y} =} gf (@var{x}, @var{m}, @var{primpoly})\n\ -+Creates a Galois field array GF(2^@var{m}) from the matrix @var{x}. The\n\ -+Galois field has 2^@var{m} elements, where @var{m} must be between 1 and 16.\n\ -+The elements of @var{x} must be between 0 and 2^@var{m} - 1. If @var{m} is\n\ -+undefined it defaults to the value 1.\n\ -+\n\ -+The primitive polynomial to use in the creation of Galois field can be\n\ -+specified with the @var{primpoly} variable. If this is undefined a default\n\ -+primitive polynomial is used. It should be noted that the primitive\n\ -+polynomial must be of the degree @var{m} and it must be irreducible.\n\ -+\n\ -+The output of this function is recognized as a Galois field by Octave and\n\ -+other matrices will be converted to the same Galois field when used in an\n\ -+arithmetic operation with a Galois field.\n\ -+\n\ -+@seealso{isprimitive, primpoly}\n\ -+@end deftypefn") -+{ -+ Matrix data; -+ octave_value retval; -+ int nargin = args.length (); -+ int m = 1; -+ int primpoly = 0; -+ -+ if (nargin < 1 || nargin > 3) -+ { -+ print_usage (); -+ return retval; -+ } -+ -+ data = args(0).matrix_value (); -+ if (nargin > 1) -+ m = args(1).int_value (); -+ if (nargin > 2) -+ primpoly = args(2).int_value (); -+ -+ if (!galois_type_loaded) -+ { -+ octave_galois::register_type (); -+ install_gm_gm_ops (); -+ install_m_gm_ops (); -+ install_gm_m_ops (); -+ install_s_gm_ops (); -+ install_gm_s_ops (); -+ galois_type_loaded = true; -+ interp.mlock (); -+ } -+ -+ retval = new octave_galois (data, m, primpoly); -+ return retval; -+} -+ -+/* -+%% Test input validation -+%!error gf () -+%!error gf (1, 2, 3, 4) -+*/ -+ -+static octave_value -+make_gdiag (const octave_value& a, const octave_value& b) -+{ -+ octave_value retval; -+ -+ if ((!galois_type_loaded) || (a.type_id () != -+ octave_galois::static_type_id ())) -+ err_wrong_type_arg ("gdiag", a); -+ else -+ { -+ galois m = ((const octave_galois&) a.get_rep ()).galois_value (); -+ int k = b.nint_value (); -+ -+ if (! error_state) -+ { -+ int nr = m.rows (); -+ int nc = m.columns (); -+ -+ if (nr == 0 || nc == 0) -+ retval = new octave_galois (m); -+ else if (nr == 1 || nc == 1) -+ { -+ int roff = 0; -+ int coff = 0; -+ if (k > 0) -+ { -+ roff = 0; -+ coff = k; -+ } -+ else if (k < 0) -+ { -+ k = -k; -+ roff = k; -+ coff = 0; -+ } -+ -+ if (nr == 1) -+ { -+ int n = nc + k; -+ galois r (n, n, 0, m.m (), m.primpoly ()); -+ for (int i = 0; i < nc; i++) -+ r (i+roff, i+coff) = m (0, i); -+ retval = new octave_galois (r); -+ } -+ else -+ { -+ int n = nr + k; -+ galois r (n, n, 0, m.m (), m.primpoly ()); -+ for (int i = 0; i < nr; i++) -+ r (i+roff, i+coff) = m (i, 0); -+ retval = new octave_galois (r); -+ } -+ } -+ else -+ { -+ galois r = m.diag (k); -+ if (r.numel () > 0) -+ retval = new octave_galois (r); -+ } -+ } -+ else -+ err_wrong_type_arg ("gdiag", a); -+ } -+ return retval; -+} -+ -+// PKG_ADD: autoload ("gdiag", "gf.oct"); -+// PKG_DEL: autoload ("gdiag", "gf.oct", "remove"); -+DEFUN_DLD (gdiag, args, , -+ "-*- texinfo -*-\n\ -+@deftypefn {Loadable Function} {} gdiag (@var{v}, @var{k})\n\ -+Return a diagonal matrix with Galois vector @var{v} on diagonal @var{k}.\n\ -+The second argument is optional. If it is positive, the vector is placed on\n\ -+the @var{k}-th super-diagonal. If it is negative, it is placed on the\n\ -+@var{-k}-th sub-diagonal. The default value of @var{k} is 0, and the\n\ -+vector is placed on the main diagonal. For example,\n\ -+\n\ -+@example\n\ -+gdiag (gf ([1, 2, 3], 2), 1)\n\ -+ans =\n\ -+GF(2^2) array. Primitive Polynomial = D^2+D+1 (decimal 7)\n\ -+\n\ -+Array elements =\n\ -+\n\ -+ 0 1 0 0\n\ -+ 0 0 2 0\n\ -+ 0 0 0 3\n\ -+ 0 0 0 0\n\ -+\n\ -+@end example\n\ -+@seealso{diag}\n\ -+@end deftypefn") -+{ -+ octave_value retval; -+ -+ int nargin = args.length (); -+ -+ if (nargin == 1 && args(0).is_defined ()) -+ retval = make_gdiag (args(0), octave_value (0.)); -+ else if (nargin == 2 && args(0).is_defined () && args(1).is_defined ()) -+ retval = make_gdiag (args(0), args(1)); -+ else -+ print_usage (); -+ -+ return retval; -+} -+ -+/* -+%% Test input validation -+%!error gdiag () -+%!error gdiag (1, 2, 3) -+*/ -+ -+// PKG_ADD: autoload ("greshape", "gf.oct"); -+// PKG_DEL: autoload ("greshape", "gf.oct", "remove"); -+DEFUN_DLD (greshape, args, , -+ "-*- texinfo -*-\n\ -+@deftypefn {Loadable Function} {} greshape (@var{a}, @var{m}, @var{n})\n\ -+Return a matrix with @var{m} rows and @var{n} columns whose elements are\n\ -+taken from the Galois array @var{a}. To decide how to order the elements,\n\ -+Octave pretends that the elements of a matrix are stored in column-major\n\ -+order (like Fortran arrays are stored).\n\ -+\n\ -+For example,\n\ -+\n\ -+@example\n\ -+greshape (gf ([1, 2, 3, 4], 3), 2, 2)\n\ -+ans =\n\ -+GF(2^3) array. Primitive Polynomial = D^3+D+1 (decimal 11)\n\ -+\n\ -+Array elements =\n\ -+\n\ -+ 1 3\n\ -+ 2 4\n\ -+\n\ -+@end example\n\ -+\n\ -+The @code{greshape} function is equivalent to\n\ -+\n\ -+@example\n\ -+@group\n\ -+retval = gf (zeros (m, n), a.m, a.prim_poly);\n\ -+retval(:) = a;\n\ -+@end group\n\ -+@end example\n\ -+\n\ -+@noindent\n\ -+but it is somewhat less cryptic to use @code{reshape} instead of the\n\ -+colon operator. Note that the total number of elements in the original\n\ -+matrix must match the total number of elements in the new matrix.\n\ -+@seealso{reshape, :}\n\ -+@end deftypefn") -+{ -+ octave_value retval; -+ int nargin = args.length (); -+ -+ if (nargin != 2 && nargin != 3) -+ { -+ print_usage (); -+ } -+ else -+ { -+ int mr = 0, mc = 0; -+ -+ if ((!galois_type_loaded) || (args(0).type_id () != -+ octave_galois::static_type_id ())) -+ { -+ err_wrong_type_arg ("greshape", args(0)); -+ return retval; -+ } -+ galois a = ((const octave_galois&) args(0).get_rep ()).galois_value (); -+ -+ if (nargin == 2) -+ { -+ RowVector tmp = args(1).row_vector_value (); -+ mr = (int)tmp(0); -+ mc = (int)tmp(1); -+ } -+ else if (nargin == 3) -+ { -+ mr = args(1).nint_value (); -+ mc = args(2).nint_value (); -+ } -+ -+ int nr = a.rows (); -+ int nc = a.cols (); -+ if ((nr * nc) != (mr * mc)) -+ error ("greshape: sizes must match"); -+ else -+ { -+ RowVector tmp1 (mr*mc); -+ for (int i = 0; i < nr; i++) -+ for (int j = 0; j < nc; j++) -+ tmp1(i+j*nr) = (double)a(i, j); -+ galois tmp2 (mr, mc, 0, a.m (), a.primpoly ()); -+ for (int i = 0; i < mr; i++) -+ for (int j = 0; j < mc; j++) -+ tmp2(i, j) = (int)tmp1(i+j*mr); -+ retval = new octave_galois (tmp2); -+ } -+ } -+ return retval; -+} -+ -+/* -+%% Test input validation -+%!error greshape () -+%!error greshape (1) -+%!error greshape (1, 2, 3, 4) -+*/ -+ -+#define DATA_REDUCTION(FCN) \ -+ \ -+ octave_value_list retval; \ -+ \ -+ int nargin = args.length (); \ -+ \ -+ if (nargin == 1 || nargin == 2) \ -+ { \ -+ octave_value arg = args(0); \ -+ \ -+ int dim = (nargin == 1 ? -1 : args(1).int_value (true) - 1); \ -+ \ -+ if (! error_state) \ -+ { \ -+ if (dim <= 1 && dim >= -1) \ -+ { \ -+ if (galois_type_loaded && (arg.type_id () == \ -+ octave_galois::static_type_id ())) \ -+ { \ -+ galois tmp = ((const octave_galois&)arg.get_rep ()).galois_value (); \ -+ \ -+ if (! error_state) \ -+ retval(0) = new octave_galois (tmp.FCN (dim)); \ -+ } \ -+ else \ -+ { \ -+ err_wrong_type_arg (#FCN, arg); \ -+ return retval; \ -+ } \ -+ } \ -+ else \ -+ error (#FCN ": invalid dimension argument = %d", dim + 1); \ -+ } \ -+ } \ -+ else \ -+ print_usage (); \ -+ \ -+ return retval -+ -+// PKG_ADD: autoload ("gprod", "gf.oct"); -+// PKG_DEL: autoload ("gprod", "gf.oct", "remove"); -+DEFUN_DLD (gprod, args, , -+ "-*- texinfo -*-\n\ -+@deftypefn {Loadable Function} {} gprod (@var{x}, @var{dim})\n\ -+Product of elements along dimension @var{dim} of Galois array. If\n\ -+@var{dim} is omitted, it defaults to 1 (column-wise products).\n\ -+@seealso{prod}\n\ -+@end deftypefn") -+{ -+ DATA_REDUCTION (prod); -+} -+ -+/* -+%% Test input validation -+%!error gprod () -+%!error gprod (1, 2, 3) -+*/ -+ -+// PKG_ADD: autoload ("gsum", "gf.oct"); -+// PKG_DEL: autoload ("gsum", "gf.oct", "remove"); -+DEFUN_DLD (gsum, args, , -+ "-*- texinfo -*-\n\ -+@deftypefn {Loadable Function} {} gsum (@var{x}, @var{dim})\n\ -+Sum of elements along dimension @var{dim} of Galois array. If @var{dim}\n\ -+is omitted, it defaults to 1 (column-wise sum).\n\ -+@seealso{sum}\n\ -+@end deftypefn") -+{ -+ DATA_REDUCTION (sum); -+} -+ -+/* -+%% Test input validation -+%!error gsum () -+%!error gsum (1, 2, 3) -+*/ -+ -+// PKG_ADD: autoload ("gsumsq", "gf.oct"); -+// PKG_DEL: autoload ("gsumsq", "gf.oct", "remove"); -+DEFUN_DLD (gsumsq, args, , -+ "-*- texinfo -*-\n\ -+@deftypefn {Loadable Function} {} gsumsq (@var{x}, @var{dim})\n\ -+Sum of squares of elements along dimension @var{dim} of Galois array.\n\ -+If @var{dim} is omitted, it defaults to 1 (column-wise sum of squares).\n\ -+\n\ -+This function is equivalent to computing\n\ -+@example\n\ -+gsum (x .* conj (x), dim)\n\ -+@end example\n\ -+but it uses less memory.\n\ -+@seealso{sumsq}\n\ -+@end deftypefn") -+{ -+ DATA_REDUCTION (sumsq); -+} -+ -+/* -+%% Test input validation -+%!error gsumsq () -+%!error gsumsq (1, 2, 3) -+*/ -+ -+// PKG_ADD: autoload ("gsqrt", "gf.oct"); -+// PKG_DEL: autoload ("gsqrt", "gf.oct", "remove"); -+DEFUN_DLD (gsqrt, args, , -+ "-*- texinfo -*-\n\ -+@deftypefn {Loadable Function} {} gsqrt (@var{x})\n\ -+Compute the square root of @var{x}, element by element, in a Galois Field.\n\ -+@seealso{exp}\n\ -+@end deftypefn") -+{ -+ octave_value retval; -+ int nargin = args.length (); -+ -+ if (nargin != 1) -+ { -+ print_usage (); -+ return retval; -+ } -+ -+ if (!galois_type_loaded || (args(0).type_id () != -+ octave_galois::static_type_id ())) -+ { -+ err_wrong_type_arg ("gsqrt", args(0)); -+ return retval; -+ } -+ -+ galois a = ((const octave_galois&) args(0).get_rep ()).galois_value (); -+ -+ retval = new octave_galois (a.sqrt ()); -+ -+ return retval; -+} -+ -+/* -+%% Test input validation -+%!error gsqrt () -+%!error gsqrt (1, 2) -+*/ -+ -+// PKG_ADD: autoload ("glog", "gf.oct"); -+// PKG_DEL: autoload ("glog", "gf.oct", "remove"); -+DEFUN_DLD (glog, args, , -+ "-*- texinfo -*-\n\ -+@deftypefn {Loadable Function} {} glog (@var{x})\n\ -+Compute the natural logarithm for each element of @var{x} for a Galois\n\ -+array.\n\ -+@seealso{log}\n\ -+@end deftypefn") -+{ -+ octave_value retval; -+ int nargin = args.length (); -+ -+ if (nargin != 1) -+ { -+ print_usage (); -+ return retval; -+ } -+ -+ if (!galois_type_loaded || (args(0).type_id () != -+ octave_galois::static_type_id ())) -+ { -+ err_wrong_type_arg ("glog", args(0)); -+ return retval; -+ } -+ -+ galois a = ((const octave_galois&) args(0).get_rep ()).galois_value (); -+ -+ retval = new octave_galois (a.log ()); -+ -+ return retval; -+} -+ -+/* -+%% Test input validation -+%!error glog () -+%!error glog (1, 2) -+*/ -+ -+// PKG_ADD: autoload ("gexp", "gf.oct"); -+// PKG_DEL: autoload ("gexp", "gf.oct", "remove"); -+DEFUN_DLD (gexp, args, , -+ "-*- texinfo -*-\n\ -+@deftypefn {Loadable Function} {} gexp (@var{x})\n\ -+Compute the anti-logarithm for each element of @var{x} for a Galois\n\ -+array.\n\ -+@seealso{exp}\n\ -+@end deftypefn") -+{ -+ octave_value retval; -+ int nargin = args.length (); -+ -+ if (nargin != 1) -+ { -+ print_usage (); -+ return retval; -+ } -+ -+ if (!galois_type_loaded || (args(0).type_id () != -+ octave_galois::static_type_id ())) -+ { -+ err_wrong_type_arg ("gexp", args(0)); -+ return retval; -+ } -+ -+ galois a = ((const octave_galois&) args(0).get_rep ()).galois_value (); -+ -+ retval = new octave_galois (a.exp ()); -+ -+ return retval; -+} -+ -+/* -+%% Test input validation -+%!error gexp () -+%!error gexp (1, 2) -+*/ -+ -+static inline int -+modn (int x, int m, int n) -+{ -+ while (x >= n) -+ { -+ x -= n; -+ x = (x >> m) + (x & n); -+ } -+ return x; -+} -+ -+galois -+filter (galois& b, galois& a, galois& x, galois& si) -+{ -+ int ab_len = (a.numel () > b.numel () ? a.numel () : b.numel ()); -+ b.resize (dim_vector (ab_len, 1), 0); -+ galois retval (x.numel (), 1, 0, b.m (), b.primpoly ()); -+ int norm = a(0, 0); -+ -+ if (norm == 0) -+ { -+ error ("gfilter: the first element of a must be non-zero"); -+ return galois (); -+ } -+ if (si.numel () != ab_len - 1) -+ { -+ error ("gfilter: si must be a vector of length max(numel(a), numel(b)) - 1"); -+ return galois (); -+ } -+ if (norm != 1) -+ { -+ int idx_norm = b.index_of (norm); -+ for (int i = 0; i < b.numel (); i++) -+ { -+ if (b(i, 0) != 0) -+ b(i, 0) = b.alpha_to (modn (b.index_of (b(i, 0))-idx_norm+b.n (), -+ b.m (), b.n ())); -+ } -+ } -+ if (a.numel () > 1) -+ { -+ a.resize (dim_vector (ab_len, 1), 0); -+ -+ if (norm != 1) -+ { -+ int idx_norm = a.index_of (norm); -+ for (int i = 0; i < a.numel (); i++) -+ if (a(i, 0) != 0) -+ a(i, 0) = a.alpha_to (modn (a.index_of (a(i, 0))-idx_norm+a.n (), -+ a.m (), a.n ())); -+ } -+ -+ for (int i = 0; i < x.numel (); i++) -+ { -+ retval(i, 0) = si(0, 0); -+ if ((b(0, 0) != 0) && (x(i, 0) != 0)) -+ retval(i, 0) ^= b.alpha_to (modn (b.index_of (b(0, 0)) + -+ b.index_of (x(i, 0)), b.m (), b.n ())); -+ if (si.numel () > 1) -+ { -+ for (int j = 0; j < si.numel () - 1; j++) -+ { -+ si(j, 0) = si(j+1, 0); -+ if ((a(j+1, 0) != 0) && (retval(i, 0) != 0)) -+ si(j, 0) ^= a.alpha_to (modn (a.index_of (a(j+1, 0)) + -+ a.index_of (retval(i, 0)), a.m (), a.n ())); -+ if ((b(j+1, 0) != 0) && (x(i, 0) != 0)) -+ si(j, 0) ^= b.alpha_to (modn (b.index_of (b(j+1, 0)) + -+ b.index_of (x(i, 0)), b.m (), b.n ())); -+ } -+ si(si.numel ()-1, 0) = 0; -+ if ((a(si.numel (), 0) != 0) && (retval(i, 0) != 0)) -+ si(si.numel ()-1, 0) ^= a.alpha_to (modn (a.index_of (a(si.numel (), 0)) -+ + a.index_of (retval(i, 0)), -+ a.m (), a.n ())); -+ if ((b(si.numel (), 0) != 0) && (x(i, 0) != 0)) -+ si(si.numel ()-1, 0) ^= b.alpha_to (modn (b.index_of (b(si.numel (), 0)) -+ + b.index_of (x(i, 0)), -+ b.m (), b.n ())); -+ } -+ else -+ { -+ si(0, 0) = 0; -+ if ((a(1, 0) != 0) && (retval(i, 0) != 0)) -+ si(0, 0) ^= a.alpha_to (modn (a.index_of (a(1, 0))+ -+ a.index_of (retval(i, 0)), a.m (), a.n ())); -+ if ((b(1, 0) != 0) && (x(i, 0) != 0)) -+ si(0, 0) ^= b.alpha_to (modn (b.index_of (b(1, 0))+ -+ b.index_of (x(i, 0)), b.m (), b.n ())); -+ } -+ } -+ } -+ else if (si.numel () > 0) -+ { -+ for (int i = 0; i < x.numel (); i++) -+ { -+ retval(i, 0) = si(0, 0); -+ if ((b(0, 0) != 0) && (x(i, 0) != 0)) -+ retval(i, 0) ^= b.alpha_to (modn (b.index_of (b(0, 0)) + -+ b.index_of (x(i, 0)), b.m (), b.n ())); -+ if (si.numel () > 1) -+ { -+ for (int j = 0; j < si.numel () - 1; j++) -+ { -+ si(j, 0) = si(j+1, 0); -+ if ((b(j+1, 0) != 0) && (x(i, 0) != 0)) -+ si(j, 0) ^= b.alpha_to (modn (b.index_of (b(j+1, 0)) + -+ b.index_of (x(i, 0)), b.m (), b.n ())); -+ } -+ si(si.numel ()-1, 0) = 0; -+ if ((b(si.numel (), 0) != 0) && (x(i, 0) != 0)) -+ si(si.numel ()-1, 0) ^= b.alpha_to (modn (b.index_of (b(si.numel (), 0)) -+ + b.index_of (x(i, 0)), -+ b.m (), b.n ())); -+ } -+ else -+ { -+ si(0, 0) = 0; -+ if ((b(1, 0) != 0) && (x(i, 0) != 0)) -+ si(0, 0) ^= b.alpha_to (modn (b.index_of (b(1, 0)) + -+ b.index_of (x(i, 0)), b.m (), b.n ())); -+ } -+ } -+ } -+ else -+ for (int i = 0; i < x.numel (); i++) -+ if ((b(0, 0) != 0) && (x(i, 0) != 0)) -+ retval(i, 0) = b.alpha_to (modn (b.index_of (b(0, 0)) + -+ b.index_of (x(i, 0)), b.m (), b.n ())); -+ -+ return retval; -+} -+ -+ -+// PKG_ADD: autoload ("gfilter", "gf.oct"); -+// PKG_DEL: autoload ("gfilter", "gf.oct", "remove"); -+DEFUN_DLD (gfilter, args, nargout, -+ "-*- texinfo -*-\n\ -+@deftypefn {Loadable Function} {y =} gfilter (@var{b}, @var{a}, @var{x})\n\ -+@deftypefnx {Loadable Function} {[@var{y}, @var{sf}] =} gfilter (@var{b}, @var{a}, @var{x}, @var{si})\n\ -+Digital filtering of vectors in a Galois Field. Returns the solution to\n\ -+the following linear, time-invariant difference equation over a Galois\n\ -+Field:\n\ -+@tex\n\ -+$$\n\ -+\\sum_{k=0}^N a_{k+1} y_{n-k} = \\sum_{k=0}^M b_{k+1} x_{n-k}, \\qquad\n\ -+ 1 \\le n \\le P\n\ -+$$\n\ -+@end tex\n\ -+@ifnottex\n\ -+\n\ -+@smallexample\n\ -+@group\n\ -+ N M\n\ -+ SUM a(k+1) y(n-k) = SUM b(k+1) x(n-k) for 1<=n<=numel(x)\n\ -+ k=0 k=0\n\ -+@end group\n\ -+@end smallexample\n\ -+@end ifnottex\n\ -+\n\ -+@noindent\n\ -+where\n\ -+@tex\n\ -+ $a \\in \\Re^{N-1}$, $b \\in \\Re^{M-1}$, and $x \\in \\Re^P$.\n\ -+@end tex\n\ -+@ifnottex\n\ -+ N=numel(a)-1 and M=numel(b)-1.\n\ -+@end ifnottex\n\ -+An equivalent form of this equation is:\n\ -+@tex\n\ -+$$\n\ -+y_n = -\\sum_{k=1}^N c_{k+1} y_{n-k} + \\sum_{k=0}^M d_{k+1} x_{n-k}, \\qquad\n\ -+ 1 \\le n \\le P\n\ -+$$\n\ -+@end tex\n\ -+@ifnottex\n\ -+\n\ -+@smallexample\n\ -+@group\n\ -+ N M\n\ -+ y(n) = - SUM c(k+1) y(n-k) + SUM d(k+1) x(n-k) for 1<=n<=numel(x)\n\ -+ k=1 k=0\n\ -+@end group\n\ -+@end smallexample\n\ -+@end ifnottex\n\ -+\n\ -+@noindent\n\ -+where\n\ -+@tex\n\ -+$c = a/a_1$ and $d = b/a_1$.\n\ -+@end tex\n\ -+@ifnottex\n\ -+ c = a/a(1) and d = b/a(1).\n\ -+@end ifnottex\n\ -+\n\ -+If the fourth argument @var{si} is provided, it is taken as the initial\n\ -+state of the system and the final state is returned as @var{sf}. The\n\ -+state vector is a column vector whose length is equal to the length of\n\ -+the longest coefficient vector minus one. If @var{si} is not supplied,\n\ -+the initial state vector is set to all zeros.\n\ -+@seealso{filter}\n\ -+@end deftypefn") -+{ -+ octave_value_list retval; -+ -+ int nargin = args.length (); -+ -+ if (nargin < 3 || nargin > 4) -+ { -+ print_usage (); -+ return retval; -+ } -+ -+ if (!galois_type_loaded) -+ { -+ error ("gfilter: wrong argument types"); -+ return retval; -+ } -+ -+ bool x_is_row_vector = (args(2).rows () == 1); -+ bool si_is_row_vector = (nargin == 4 && args(3).rows () == 1); -+ galois b, a, x, si; -+ bool ib=false, ia=false, ix = false, isi=false; -+ -+ if (args(0).type_id () == octave_galois::static_type_id ()) -+ { -+ b = ((const octave_galois&) args(0).get_rep ()).galois_value (); -+ ib = true; -+ } -+ if (args(1).type_id () == octave_galois::static_type_id ()) -+ { -+ a = ((const octave_galois&) args(1).get_rep ()).galois_value (); -+ ia = true; -+ } -+ if (args(2).type_id () == octave_galois::static_type_id ()) -+ { -+ x = ((const octave_galois&) args(2).get_rep ()).galois_value (); -+ ix = true; -+ } -+ if (nargin == 4) -+ { -+ if (args(3).type_id () == octave_galois::static_type_id ()) -+ { -+ si = ((const octave_galois&) args(3).get_rep ()).galois_value (); -+ isi = true; -+ } -+ } -+ -+ if (!ib && !ia && !ix && !isi) -+ { -+ error ("gfilter: wrong argument types"); -+ return retval; -+ } -+ -+ if (!ib) -+ { -+ if (ia) -+ b = galois (args(0).matrix_value (), a.m (), a.primpoly ()); -+ else if (ix) -+ b = galois (args(0).matrix_value (), x.m (), x.primpoly ()); -+ else if (isi) -+ b = galois (args(0).matrix_value (), si.m (), si.primpoly ()); -+ } -+ if (!ia) -+ a = galois (args(1).matrix_value (), b.m (), b.primpoly ()); -+ if (!ix) -+ x = galois (args(2).matrix_value (), b.m (), b.primpoly ()); -+ -+ if (nargin == 4) -+ { -+ if (!isi) -+ si = galois (args(3).matrix_value (), b.m (), b.primpoly ()); -+ } -+ else -+ { -+ int a_len = a.numel (); -+ int b_len = b.numel (); -+ -+ int si_len = (a_len > b_len ? a_len : b_len) - 1; -+ -+ si = galois (si_len, 1, 0, b.m (), b.primpoly ()); -+ } -+ -+ if ((b.m () != a.m ()) || (b.m () != x.m ()) || (b.m () != si.m ()) || -+ (b.primpoly () != a.primpoly ()) || (b.primpoly () != x.primpoly ()) || -+ (b.primpoly () != si.primpoly ())) -+ { -+ error ("gfilter: arguments must be in same galois field"); -+ return retval; -+ } -+ -+ if (b.cols () > 1) -+ b = b.transpose (); -+ if (a.cols () > 1) -+ a = a.transpose (); -+ if (x.cols () > 1) -+ x = x.transpose (); -+ if (si.cols () > 1) -+ si = si.transpose (); -+ -+ if (b.cols () > 1 || a.cols () > 1 || x.cols () > 1 || si.cols () > 1) -+ { -+ error ("gfilter: arguments must be vectors"); -+ return retval; -+ } -+ -+ galois y (filter (b, a, x, si)); -+ if (nargout == 2) -+ { -+ if (si_is_row_vector) -+ retval(1) = new octave_galois (si.transpose ()); -+ else -+ retval(1) = new octave_galois (si); -+ } -+ -+ if (x_is_row_vector) -+ retval(0) = new octave_galois (y.transpose ()); -+ else -+ retval(0) = new octave_galois (y); -+ -+ return retval; -+} -+ -+/* -+%% Test input validation -+%!error gfilter () -+%!error gfilter (1) -+%!error gfilter (1, 2) -+%!error gfilter (1, 2, 3, 4, 5) -+*/ -+ -+// PKG_ADD: autoload ("glu", "gf.oct"); -+// PKG_DEL: autoload ("glu", "gf.oct", "remove"); -+DEFUN_DLD (glu, args, nargout, -+ "-*- texinfo -*-\n\ -+@deftypefn {Loadable Function} {[@var{l}, @var{u}, @var{p}] =} glu (@var{a})\n\ -+@cindex LU decomposition of Galois matrix\n\ -+Compute the LU decomposition of @var{a} in a Galois Field. The result is\n\ -+returned in a permuted form, according to the optional return value\n\ -+@var{p}. For example, given the matrix\n\ -+@code{a = gf ([1, 2; 3, 4], 3)},\n\ -+\n\ -+@example\n\ -+[l, u, p] = glu (a)\n\ -+@end example\n\ -+\n\ -+@noindent\n\ -+returns\n\ -+\n\ -+@example\n\ -+l =\n\ -+GF(2^3) array. Primitive Polynomial = D^3+D+1 (decimal 11)\n\ -+\n\ -+Array elements =\n\ -+\n\ -+ 1 0\n\ -+ 6 1\n\ -+\n\ -+u =\n\ -+GF(2^3) array. Primitive Polynomial = D^3+D+1 (decimal 11)\n\ -+\n\ -+Array elements =\n\ -+\n\ -+ 3 4\n\ -+ 0 7\n\ -+\n\ -+p =\n\ -+\n\ -+Permutation Matrix\n\ -+\n\ -+ 0 1\n\ -+ 1 0\n\ -+\n\ -+@end example\n\ -+\n\ -+Such that @code{@var{p} * @var{a} = @var{l} * @var{u}}. If the argument\n\ -+@var{p} is not included then the permutations are applied to @var{l}\n\ -+so that @code{@var{a} = @var{l} * @var{u}}. @var{l} is then a pseudo-\n\ -+lower triangular matrix. The matrix @var{a} can be rectangular.\n\ -+@seealso{lu}\n\ -+@end deftypefn") -+{ -+ octave_value_list retval; -+ -+ -+ int nargin = args.length (); -+ -+ if (nargin != 1 || nargout > 3) -+ { -+ print_usage (); -+ return retval; -+ } -+ -+ octave_value arg = args(0); -+ -+ if (!galois_type_loaded || (arg.type_id () != -+ octave_galois::static_type_id ())) -+ { -+ err_wrong_type_arg ("glu", arg); -+ return retval; -+ } -+ -+ galois m = ((const octave_galois&) arg.get_rep ()).galois_value (); -+ -+ int nr = arg.rows (); -+ int nc = arg.columns (); -+ -+ int arg_is_empty = arg.isempty (); -+ -+ if (arg_is_empty < 0) -+ return retval; -+ else if (arg_is_empty > 0) -+ { -+ retval(0) = new octave_galois (galois (0, 0, 0, m.m (), m.primpoly ())); -+ retval(1) = new octave_galois (galois (0, 0, 0, m.m (), m.primpoly ())); -+ retval(2) = new octave_galois (galois (0, 0, 0, m.m (), m.primpoly ())); -+ return retval; -+ } -+ -+ if (! error_state) -+ { -+ galoisLU fact (m); -+ -+ switch (nargout) -+ { -+ case 0: -+ case 1: -+ case 2: -+ { -+ // While we don't have sparse galois matrices converting the -+ // permutation matrix to a full matrix is the best we can do. -+ Matrix P = Matrix (fact.P ()); -+ galois L = P.transpose () * fact.L (); -+ retval(1) = new octave_galois (fact.U ()); -+ retval(0) = new octave_galois (L); -+ } -+ break; -+ -+ case 3: -+ default: -+ retval(2) = fact.P (); -+ retval(1) = new octave_galois (fact.U ()); -+ retval(0) = new octave_galois (fact.L ()); -+ break; -+ } -+ } -+ -+ return retval; -+} -+ -+/* -+%% Test input validation -+%!error glu () -+%!error glu (1, 2) -+*/ -+ -+// PKG_ADD: autoload ("ginv", "gf.oct"); -+// PKG_DEL: autoload ("ginv", "gf.oct", "remove"); -+DEFUN_DLD (ginv, args, nargout, -+ "-*- texinfo -*-\n\ -+@deftypefn {Loadable Function} {[@var{x}, @var{rcond}] =} ginv (@var{a})\n\ -+Compute the inverse of the square matrix @var{a}. Return an estimate\n\ -+of the reciprocal condition number if requested, otherwise warn of an\n\ -+ill-conditioned matrix if the reciprocal condition number is small.\n\ -+@seealso{inv}\n\ -+@end deftypefn") -+{ -+ octave_value_list retval; -+ -+ int nargin = args.length (); -+ -+ if (nargin != 1) -+ { -+ print_usage (); -+ return retval; -+ } -+ -+ octave_value arg = args(0); -+ -+ int nr = arg.rows (); -+ int nc = arg.columns (); -+ -+ if (!galois_type_loaded || (arg.type_id () != -+ octave_galois::static_type_id ())) -+ { -+ err_wrong_type_arg ("ginverse", arg); -+ return retval; -+ } -+ -+ galois m = ((const octave_galois&) arg.get_rep ()).galois_value (); -+ -+ int arg_is_empty = arg.isempty (); -+ -+ if (arg_is_empty < 0) -+ return retval; -+ else if (arg_is_empty > 0) -+ { -+ retval(0) = new octave_galois (galois (0, 0, 0, m.m (), m.primpoly ())); -+ return retval; -+ } -+ if (nr != nc) -+ { -+ err_square_matrix_required ("ginverse", "X"); -+ return retval; -+ } -+ -+ if (! error_state) -+ { -+ int info; -+ double rcond = 0.0; -+ -+ galois result = m.inverse (info, 1); -+ -+ if (nargout > 1) -+ retval(1) = rcond; -+ -+ retval(0) = new octave_galois (result); -+ -+ if (nargout < 2 && info == -1) -+ warning ("inverse: matrix singular to machine precision, rcond = %g", rcond); -+ } -+ -+ return retval; -+} -+ -+/* -+%% Test input validation -+%!error ginv () -+%!error ginv (1, 2) -+*/ -+ -+// FIXME: this should really be done with an alias, but -+// alias_builtin() won't do the right thing if we are actually using -+// dynamic linking. -+ -+// PKG_ADD: autoload ("ginverse", "gf.oct"); -+// PKG_DEL: autoload ("ginverse", "gf.oct", "remove"); -+DEFUN_DLD (ginverse, args, nargout, -+ "-*- texinfo -*-\n\ -+@deftypefn {Loadable Function} {} ginverse (@var{a})\n\ -+Compute the inverse of the square matrix @var{a}. Return an estimate\n\ -+of the reciprocal condition number if requested, otherwise warn of an\n\ -+ill-conditioned matrix if the reciprocal condition number is small.\n\ -+@seealso{ginv}\n\ -+@end deftypefn") -+{ -+ return Fginv (args, nargout); -+} -+ -+/* -+%% Test input validation -+%!error ginverse () -+%!error ginverse (1, 2) -+*/ -+ -+// PKG_ADD: autoload ("gdet", "gf.oct"); -+// PKG_DEL: autoload ("gdet", "gf.oct", "remove"); -+DEFUN_DLD (gdet, args, nargout, -+ "-*- texinfo -*-\n\ -+@deftypefn {Loadable Function} {@var{d} =} gdet (@var{a})\n\ -+Compute the determinant of the Galois array @var{a}.\n\ -+@seealso{det}\n\ -+@end deftypefn") -+{ -+ octave_value retval; -+ -+ int nargin = args.length (); -+ -+ if (nargin != 1) -+ { -+ print_usage (); -+ return retval; -+ } -+ -+ octave_value arg = args(0); -+ -+ if (!galois_type_loaded || (arg.type_id () != -+ octave_galois::static_type_id ())) -+ { -+ err_wrong_type_arg ("gdet", arg); -+ return retval; -+ } -+ -+ int nr = arg.rows (); -+ int nc = arg.columns (); -+ -+ galois m = ((const octave_galois&) arg.get_rep ()).galois_value (); -+ -+ int arg_is_empty = arg.isempty (); -+ -+ if (arg_is_empty < 0) -+ return retval; -+ else if (arg_is_empty > 0) -+ { -+ retval = new octave_galois (galois (1, 1, 1, m.m (), m.primpoly ())); -+ return retval; -+ } -+ -+ if (nr != nc) -+ { -+ err_square_matrix_required ("det", "A"); -+ return retval; -+ } -+ -+ retval = new octave_galois (m.determinant ()); -+ return retval; -+} -+ -+/* -+%% Test input validation -+%!error gdet () -+%!error gdet (1, 2) -+*/ -+ -+// PKG_ADD: autoload ("grank", "gf.oct"); -+// PKG_DEL: autoload ("grank", "gf.oct", "remove"); -+DEFUN_DLD (grank, args, nargout, -+ "-*- texinfo -*-\n\ -+@deftypefn {Loadable Function} {@var{d} =} grank (@var{a})\n\ -+Compute the rank of the Galois array @var{a} by counting the independent\n\ -+rows and columns.\n\ -+@seealso{rank}\n\ -+@end deftypefn") -+{ -+ octave_value retval; -+ -+ int nargin = args.length (); -+ -+ if (nargin != 1) -+ { -+ print_usage (); -+ return retval; -+ } -+ -+ octave_value arg = args(0); -+ -+ if (!galois_type_loaded || (arg.type_id () != -+ octave_galois::static_type_id ())) -+ { -+ err_wrong_type_arg ("grank", arg); -+ return retval; -+ } -+ -+ int nr = arg.rows (); -+ int nc = arg.columns (); -+ -+ galois m = ((const octave_galois&) arg.get_rep ()).galois_value (); -+ -+ int arg_is_empty = arg.isempty (); -+ -+ if (arg_is_empty > 0) -+ retval = 0.0; -+ else if (arg_is_empty == 0) -+ { -+ int d = 0; -+ int mm = m.m (); -+ int mn = m.n (); -+ OCTAVE_LOCAL_BUFFER (int, ci, nr); -+ -+ for (int i = 0; i < nc; i++) -+ { -+ int idx = -1; -+ int iel = 0; -+ for (int j = 0; j < nr; j++) -+ { -+ ci[j] = m.elem (j, i); -+ if (ci[j] != 0 && idx == -1) -+ { -+ iel = ci[j]; -+ idx = j; -+ } -+ } -+ -+ if (idx != -1) -+ { -+ d++; -+ int indx = m.index_of (iel); -+ for (int j = 0; j < nr; j++) -+ if (ci[j] != 0) -+ ci[j] = m.alpha_to (modn (m.index_of (ci[j]) - indx + mn, mm, mn)); -+ -+ for (int j = i+1; j < nc; j++) -+ { -+ if (m.elem (idx, j) != 0) -+ { -+ indx = m.index_of (m.elem (idx, j)); -+ for (int k = 0; k < nr; k++) -+ if (ci[k] != 0) -+ m.elem (k, j) ^= m.alpha_to (modn (m.index_of (ci[k]) + indx + -+ mn, mm, mn)); -+ } -+ } -+ } -+ } -+ retval = (double)d; -+ } -+ return retval; -+} -+ -+/* -+%% Test input validation -+%!error grank () -+%!error grank (1, 2) -+*/ -+ -+// PKG_ADD: autoload ("rsenc", "gf.oct"); -+// PKG_DEL: autoload ("rsenc", "gf.oct", "remove"); -+DEFUN_DLD (rsenc, args, nargout, -+ "-*- texinfo -*-\n\ -+@deftypefn {Loadable Function} {@var{code} =} rsenc (@var{msg}, @var{n}, @var{k})\n\ -+@deftypefnx {Loadable Function} {@var{code} =} rsenc (@var{msg}, @var{n}, @var{k}, @var{g})\n\ -+@deftypefnx {Loadable Function} {@var{code} =} rsenc (@var{msg}, @var{n}, @var{k}, @var{fcr}, @var{prim})\n\ -+@deftypefnx {Loadable Function} {@var{code} =} rsenc (@dots{}, @var{parpos})\n\ -+Encodes the message @var{msg} using a [@var{n},@var{k}] Reed-Solomon coding.\n\ -+The variable @var{msg} is a Galois array with @var{k} columns and an arbitrary\n\ -+number of rows. Each row of @var{msg} represents a single block to be coded\n\ -+by the Reed-Solomon coder. The coded message is returned in the Galois\n\ -+array @var{code} containing @var{n} columns and the same number of rows as\n\ -+@var{msg}.\n\ -+\n\ -+The use of @code{rsenc} can be seen in the following short example.\n\ -+\n\ -+@example\n\ -+m = 3; n = 2^m -1; k = 3;\n\ -+msg = gf ([1 2 3; 4 5 6], m);\n\ -+code = rsenc (msg, n, k);\n\ -+@end example\n\ -+\n\ -+If @var{n} does not equal @code{2^@var{m}-1}, where m is an integer, then a\n\ -+shorten Reed-Solomon coding is used where zeros are added to the start of\n\ -+each row to obtain an allowable codeword length. The returned @var{code}\n\ -+has these prepending zeros stripped.\n\ -+\n\ -+By default the generator polynomial used in the Reed-Solomon coding is based\n\ -+on the properties of the Galois Field in which @var{msg} is given. This\n\ -+default generator polynomial can be overridden by a polynomial in @var{g}.\n\ -+Suitable generator polynomials can be constructed with @code{rsgenpoly}.\n\ -+@var{fcr} is an integer value, and it is taken to be the first consecutive\n\ -+root of the generator polynomial. The variable @var{prim} is then the\n\ -+primitive element used to construct the generator polynomial, such that\n\ -+@tex\n\ -+$g = (x - A^b) (x - A^{b+p}) \\cdots (x - A ^{b+2tp-1})$.\n\ -+@end tex\n\ -+@ifnottex\n\ -+\n\ -+@var{g} = (@var{x} - A^@var{b}) * (@var{x} - A^(@var{b}+@var{prim})) * ... * (@var{x} - A^(@var{b}+2*@var{t}*@var{prim}-1)).\n\ -+@end ifnottex\n\ -+\n\ -+where @var{b} is equal to @code{@var{fcr} * @var{prim}}. By default @var{fcr}\n\ -+and @var{prim} are both 1.\n\ -+\n\ -+By default the parity symbols are placed at the end of the coded message.\n\ -+The variable @var{parpos} controls this positioning and can take the values\n\ -+@code{\"beginning\"} or @code{\"end\"}.\n\ -+@seealso{gf, rsdec, rsgenpoly}\n\ -+@end deftypefn") -+{ -+ octave_value retval; -+ int nargin = args.length (); -+ -+ if (nargin < 3 || nargin > 5) -+ { -+ print_usage (); -+ return retval; -+ } -+ -+ if (!galois_type_loaded || (args(0).type_id () != -+ octave_galois::static_type_id ())) -+ { -+ err_wrong_type_arg ("rsenc", args(0)); -+ return retval; -+ } -+ -+ galois msg = ((const octave_galois&) args(0).get_rep ()).galois_value (); -+ int nsym = msg.rows (); -+ int primpoly = msg.primpoly (); -+ int n = args(1).nint_value (); -+ int k = args(2).nint_value (); -+ -+ int m = 1; -+ while (n > (1< __OCTAVE_GALOIS_MAX_M)) -+ { -+ error ("rsenc: invalid values of message and codeword length"); -+ return retval; -+ } -+ -+ if ((n-k) & 1) -+ { -+ error ("rsenc: difference of message and codeword length must be even"); -+ return retval; -+ } -+ -+ int nroots = n-k; -+ galois genpoly; -+ bool have_genpoly = false; -+ bool parity_at_end = true; -+ int fcr = 0; -+ int prim = 0; -+ -+ for (int i = 3; i < nargin; i++) -+ { -+ if (args(i).is_string ()) -+ { -+ std::string parstr = args(i).string_value (); -+ for (int j = 0; j < (int)parstr.length (); j++) -+ parstr[j] = toupper (parstr[j]); -+ -+ if (!parstr.compare("END")) -+ { -+ parity_at_end = true; -+ } -+ else if (!parstr.compare("BEGINNING")) -+ { -+ parity_at_end = false; -+ } -+ else -+ { -+ error ("rsenc: unrecoginized parity position"); -+ return retval; -+ } -+ } -+ else -+ { -+ if (args(i).type_id () == octave_galois::static_type_id ()) -+ { -+ if (have_genpoly) -+ { -+ print_usage (); -+ return retval; -+ } -+ genpoly = ((const octave_galois&) args(i).get_rep ()).galois_value (); -+ -+ if (genpoly.cols () > genpoly.rows ()) -+ genpoly = genpoly.transpose (); -+ } -+ else -+ { -+ if (have_genpoly) -+ { -+ if (prim != 0) -+ { -+ print_usage (); -+ return retval; -+ } -+ prim = args(i).nint_value (); -+ } -+ else -+ fcr = args(i).nint_value (); -+ } -+ have_genpoly = true; -+ } -+ } -+ -+ if ((genpoly.rows () == 0) || (genpoly.cols () == 0)) -+ { -+ if (fcr == 0) -+ fcr = 1; -+ if (prim == 0) -+ prim = 1; -+ -+ // Create polynomial of right length. -+ genpoly = galois (nroots+1, 1, 0, m, primpoly); -+ -+ genpoly(nroots, 0) = 1; -+ int i, root; -+ for (i = 0, root=fcr*prim; i < nroots; i++, root += prim) -+ { -+ genpoly(nroots-i-1, 0) = 1; -+ -+ // Multiply genpoly by @**(root + x) -+ for (int j = i; j > 0; j--) -+ { -+ int k = nroots - j; -+ if (genpoly(k, 0) != 0) -+ genpoly(k, 0) = genpoly(k+1, 0) -+ ^ genpoly.alpha_to (modn (genpoly.index_of (genpoly(k, 0)) -+ + root, m, n)); -+ else -+ genpoly(k, 0) = genpoly(k+1, 0); -+ } -+ // genpoly(nroots,0) can never be zero -+ genpoly(nroots, 0) = genpoly.alpha_to (modn (genpoly.index_of (genpoly(nroots, 0)) -+ + root, m, n)); -+ } -+ -+ } -+ else -+ { -+ if (genpoly.cols () != 1) -+ { -+ error ("rsenc: the generator polynomial must be a vector"); -+ return retval; -+ } -+ -+ if (genpoly.primpoly () != primpoly) -+ { -+ error ("rsenc: the generator polynomial must be same galois field " -+ "as the message"); -+ return retval; -+ } -+ -+ if (genpoly.rows () != nroots+1) -+ { -+ error ("rsenc: generator polynomial has incorrect order"); -+ return retval; -+ } -+ } -+ -+ int norm = genpoly(0, 0); -+ -+ // Take logarithm of generator polynomial, for faster coding -+ for (int i = 0; i < nroots+1; i++) -+ genpoly(i, 0) = genpoly.index_of (genpoly(i, 0)); -+ -+ // Add space for parity block -+ msg.resize (dim_vector (nsym, n), 0); -+ -+ // The code below basically finds the parity bits by treating the -+ // message as a polynomial and dividing it by the generator polynomial. -+ // The parity bits are then the remainder of this division. If the parity -+ // is at the end the polynomial is treat MSB first, otherwise it is -+ // treated LSB first -+ // -+ // This code could just as easily be written as -+ // [ignore par] = gdeconv(msg, genpoly); -+ // But the code below has the advantage of being 20 times faster :-) -+ -+ if (parity_at_end) -+ { -+ for (int l = 0; l < nsym; l++) -+ { -+ galois par (nroots, 1, 0, m, primpoly); -+ for (int i = 0; i < k; i++) -+ { -+ int feedback = par.index_of (par(0, 0) ^ msg(l, i)); -+ if (feedback != nn) -+ { -+ if (norm != 1) -+ feedback = modn (nn-genpoly(0, 0)+feedback, m, nn); -+ for (int j = 1; j < nroots; j++) -+ par(j, 0) ^= par.alpha_to (modn (feedback + -+ genpoly(j, 0), m, nn)); -+ } -+ for (int j = 1; j < nroots; j++) -+ par(j-1, 0) = par(j, 0); -+ if (feedback != nn) -+ par(nroots-1, 0) = par.alpha_to (modn (feedback+ -+ genpoly(nroots, 0), m, nn)); -+ else -+ par(nroots-1, 0) = 0; -+ } -+ for (int j = 0; j < nroots; j++) -+ msg(l, k+j) = par(j, 0); -+ } -+ } -+ else -+ { -+ for (int l = 0; l < nsym; l++) -+ { -+ for (int i=k; i > 0; i--) -+ msg(l, i+nroots-1) = msg(l, i-1); -+ for (int i = 0; i nroots; i--) -+ { -+ int feedback = par.index_of (par(0, 0) ^ msg(l, i-1)); -+ if (feedback != nn) -+ { -+ if (norm != 1) -+ feedback = modn (nn-genpoly(0, 0)+feedback, m, nn); -+ for (int j = 1; j < nroots; j++) -+ par(j, 0) ^= par.alpha_to (modn (feedback + -+ genpoly(j, 0), m, nn)); -+ } -+ for (int j = 1; j < nroots; j++) -+ par(j-1, 0) = par(j, 0); -+ if (feedback != nn) -+ par(nroots-1, 0) = par.alpha_to (modn (feedback+ -+ genpoly(nroots, 0), m, nn)); -+ else -+ par(nroots-1, 0) = 0; -+ } -+ for (int j = 0; j < nroots; j++) -+ msg(l, j) = par(nroots-j-1, 0); -+ } -+ } -+ -+ retval = new octave_galois (msg); -+ -+ return retval; -+} -+ -+/* -+%% Test input validation -+%!error rsenc () -+%!error rsenc (1) -+%!error rsenc (1, 2) -+%!error rsenc (1, 2, 3, 4, 5, 6) -+*/ -+ -+int -+decode_rs(galois& data, const int prim, const int iprim, const int nroots, -+ const int fcr, const int drow, const bool msb_first) -+{ -+ int deg_lambda, el, deg_omega; -+ int i, j, r, k; -+ int q, tmp, num1, num2, den, discr_r; -+ int syn_error, count; -+ int m = data.m (); -+ int n = data.n (); -+ int A0 = n; -+ -+ /* Err Locator and syndrome poly */ -+ OCTAVE_LOCAL_BUFFER (int, lambda, nroots+1); -+ OCTAVE_LOCAL_BUFFER (int, s, nroots); -+ -+ OCTAVE_LOCAL_BUFFER (int, b, nroots+1); -+ OCTAVE_LOCAL_BUFFER (int, t, nroots+1); -+ OCTAVE_LOCAL_BUFFER (int, omega, nroots+1); -+ -+ OCTAVE_LOCAL_BUFFER (int, root, nroots); -+ OCTAVE_LOCAL_BUFFER (int, reg, nroots+1); -+ OCTAVE_LOCAL_BUFFER (int, loc, nroots); -+ -+ /* form the syndromes; i.e., evaluate data(x) at roots of g(x) */ -+ if (msb_first) -+ { -+ for (i = 0; i < nroots; i++) -+ s[i] = data(drow, 0); -+ -+ for (j = 1; j < n; j++) -+ for (i = 0; i0; j--) -+ for (i = 0; i < nroots; i++) -+ if(s[i] == 0) -+ s[i] = data(drow, j-1); -+ else -+ s[i] = data(drow, j-1) ^ data.alpha_to (modn (data.index_of (s[i]) + -+ (fcr+i)*prim, m, n)); -+ } -+ -+ /* Convert syndromes to index form, checking for nonzero condition */ -+ syn_error = 0; -+ for (i = 0; i < nroots; i++) -+ { -+ syn_error |= s[i]; -+ s[i] = data.index_of (s[i]); -+ } -+ -+ if (!syn_error) -+ /* if syndrome is zero, data(drow,:) is a codeword and there are no -+ * errors to correct. So return data(drow,:) unmodified -+ */ -+ return 0; -+ -+ memset(&lambda[1], 0, nroots*sizeof (lambda[0])); -+ lambda[0] = 1; -+ -+ for (i = 0; i < nroots+1; i++) -+ b[i] = data.index_of (lambda[i]); -+ -+ /* -+ * Begin Berlekamp-Massey algorithm to determine error locator polynomial -+ */ -+ r = 0; -+ el = 0; -+ while (++r <= nroots) -+ {/* r is the step number */ -+ /* Compute discrepancy at the r-th step in poly-form */ -+ discr_r = 0; -+ for (i = 0; i < r; i++) -+ { -+ if ((lambda[i] != 0) && (s[r-i-1] != A0)) -+ { -+ discr_r ^= data.alpha_to (modn (data.index_of (lambda[i]) + -+ s[r-i-1], m, n)); -+ } -+ } -+ discr_r = data.index_of (discr_r); /* Index form */ -+ if (discr_r == A0) -+ { -+ /* 2 lines below: B(x) <-- x*B(x) */ -+ memmove(&b[1], b, nroots*sizeof (b[0])); -+ b[0] = A0; -+ } -+ else -+ { -+ /* 7 lines below: T(x) <-- lambda(x) - discr_r*x*b(x) */ -+ t[0] = lambda[0]; -+ for (i = 0 ; i < nroots; i++) -+ { -+ if(b[i] != A0) -+ t[i+1] = lambda[i+1] ^ data.alpha_to (modn (discr_r + b[i], m, n)); -+ else -+ t[i+1] = lambda[i+1]; -+ } -+ if (2 * el <= r - 1) -+ { -+ el = r - el; -+ /* -+ * 2 lines below: B(x) <-- inv(discr_r) * -+ * lambda(x) -+ */ -+ for (i = 0; i <= nroots; i++) -+ b[i] = (lambda[i] == 0) ? A0 : modn (data.index_of (lambda[i]) - -+ discr_r + n, m, n); -+ } -+ else -+ { -+ /* 2 lines below: B(x) <-- x*B(x) */ -+ memmove(&b[1], b, nroots*sizeof (b[0])); -+ b[0] = A0; -+ } -+ memcpy(lambda, t, (nroots+1)*sizeof (t[0])); -+ } -+ } -+ -+ /* Convert lambda to index form and compute deg(lambda(x)) */ -+ deg_lambda = 0; -+ for (i = 0; i < nroots+1; i++) -+ { -+ lambda[i] = data.index_of (lambda[i]); -+ if(lambda[i] != A0) -+ deg_lambda = i; -+ } -+ -+ /* Find roots of the error locator polynomial by Chien search */ -+ memcpy(®[1], &lambda[1], nroots*sizeof (reg[0])); -+ count = 0; /* Number of roots of lambda(x) */ -+ for (i = 1, k = iprim-1; i <= n; i++, k = modn (k+iprim, m, n)) -+ { -+ q = 1; /* lambda[0] is always 0 */ -+ for (j = deg_lambda; j > 0; j--) -+ { -+ if (reg[j] != A0) -+ { -+ reg[j] = modn (reg[j] + j, m, n); -+ q ^= data.alpha_to (reg[j]); -+ } -+ } -+ if (q != 0) -+ continue; /* Not a root */ -+ /* store root (index-form) and error location number */ -+ root[count] = i; -+ loc[count] = k; -+ /* If we've already found max possible roots, -+ * abort the search to save time -+ */ -+ if(++count == deg_lambda) -+ break; -+ } -+ if (deg_lambda != count) -+ { -+ /* -+ * deg(lambda) unequal to number of roots => uncorrectable -+ * error detected -+ */ -+ return -1; -+ } -+ /* -+ * Compute err evaluator poly omega(x) = s(x)*lambda(x) (modulo -+ * x**nroots). in index form. Also find deg(omega). -+ */ -+ deg_omega = 0; -+ for (i = 0; i < nroots; i++) -+ { -+ tmp = 0; -+ j = (deg_lambda < i) ? deg_lambda : i; -+ for (; j >= 0; j--) -+ { -+ if ((s[i - j] != A0) && (lambda[j] != A0)) -+ tmp ^= data.alpha_to (modn (s[i - j] + lambda[j], m, n)); -+ } -+ if(tmp != 0) -+ deg_omega = i; -+ omega[i] = data.index_of (tmp); -+ } -+ omega[nroots] = A0; -+ -+ /* -+ * Compute error values in poly-form. num1 = omega(inv(X(l))), num2 = -+ * inv(X(l))**(fcr-1) and den = lambda_pr(inv(X(l))) all in poly-form -+ */ -+ for (j = count-1; j >= 0; j--) -+ { -+ num1 = 0; -+ for (i = deg_omega; i >= 0; i--) -+ { -+ if (omega[i] != A0) -+ num1 ^= data.alpha_to (modn (omega[i] + i * root[j], m, n)); -+ } -+ num2 = data.alpha_to (modn (root[j] * (fcr - 1) + n, m, n)); -+ den = 0; -+ -+ /* lambda[i+1] for i even is the formal deriv lambda_pr of lambda[i] */ -+ for (i = (deg_lambda < nroots-1 ? deg_lambda : nroots-1) & ~1; i >= 0; -+ i -=2) -+ { -+ if(lambda[i+1] != A0) -+ den ^= data.alpha_to (modn (lambda[i+1] + i * root[j], m, n)); -+ } -+ if (den == 0) -+ { -+ count = -1; -+ break; -+ } -+ /* Apply error to data */ -+ if (num1 != 0) -+ { -+ if (msb_first) -+ data(drow, loc[j]) ^= data.alpha_to (modn (data.index_of (num1) -+ + data.index_of (num2) -+ + n - data.index_of (den), -+ m, n)); -+ else -+ data(drow, n-loc[j]-1) ^= data.alpha_to (modn (data.index_of (num1) -+ + data.index_of (num2) -+ + n - data.index_of (den), -+ m, n)); -+ } -+ } -+ -+ return count; -+} -+ -+// PKG_ADD: autoload ("rsdec", "gf.oct"); -+// PKG_DEL: autoload ("rsdec", "gf.oct", "remove"); -+DEFUN_DLD (rsdec, args, nargout, -+ "-*- texinfo -*-\n\ -+@deftypefn {Loadable Function} {@var{msg} =} rsdec (@var{code}, @var{n}, @var{k})\n\ -+@deftypefnx {Loadable Function} {@var{msg} =} rsdec (@var{code}, @var{n}, @var{k}, @var{g})\n\ -+@deftypefnx {Loadable Function} {@var{msg} =} rsdec (@var{code}, @var{n}, @var{k}, @var{fcr}, @var{prim})\n\ -+@deftypefnx {Loadable Function} {@var{msg} =} rsdec (@dots{}, @var{parpos})\n\ -+@deftypefnx {Loadable Function} {[@var{msg}, @var{nerr}] =} rsdec (@dots{})\n\ -+@deftypefnx {Loadable Function} {[@var{msg}, @var{nerr}, @var{ccode}] =} rsdec (@dots{})\n\ -+Decodes the message contained in @var{code} using a [@var{n},@var{k}]\n\ -+Reed-Solomon code. The variable @var{code} must be a Galois array with\n\ -+@var{n} columns and an arbitrary number of rows. Each row of @var{code}\n\ -+represents a single block to be decoded by the Reed-Solomon coder. The\n\ -+decoded message is returned in the variable @var{msg} containing @var{k}\n\ -+columns and the same number of rows as @var{code}.\n\ -+\n\ -+If @var{n} does not equal @code{2^@var{m}-1}, where m is an integer, then a\n\ -+shorten Reed-Solomon decoding is used where zeros are added to the start of\n\ -+each row to obtain an allowable codeword length. The returned @var{msg}\n\ -+has these prepending zeros stripped.\n\ -+\n\ -+By default the generator polynomial used in the Reed-Solomon coding is based\n\ -+on the properties of the Galois Field in which @var{msg} is given. This\n\ -+default generator polynomial can be overridden by a polynomial in @var{g}.\n\ -+Suitable generator polynomials can be constructed with @code{rsgenpoly}.\n\ -+@var{fcr} is an integer value, and it is taken to be the first consecutive\n\ -+root of the generator polynomial. The variable @var{prim} is then the\n\ -+primitive element used to construct the generator polynomial. By default\n\ -+@var{fcr} and @var{prim} are both 1. It is significantly faster to specify\n\ -+the generator polynomial in terms of @var{fcr} and @var{prim}, since @var{g}\n\ -+is converted to this form in any case.\n\ -+\n\ -+By default the parity symbols are placed at the end of the coded message.\n\ -+The variable @var{parpos} controls this positioning and can take the values\n\ -+@code{\"beginning\"} or @code{\"end\"}. If the parity symbols are at the end, the message is\n\ -+treated with the most-significant symbol first, otherwise the message is\n\ -+treated with the least-significant symbol first.\n\ -+@seealso{gf, rsenc, rsgenpoly}\n\ -+@end deftypefn") -+{ -+ octave_value_list retval; -+ -+ int nargin = args.length (); -+ -+ if (nargin < 3 || nargin > 5) -+ { -+ print_usage (); -+ return retval; -+ } -+ -+ if (!galois_type_loaded || (args(0).type_id () != -+ octave_galois::static_type_id ())) -+ { -+ err_wrong_type_arg ("rsdec", args(0)); -+ return retval; -+ } -+ -+ galois code = ((const octave_galois&) args(0).get_rep ()).galois_value (); -+ int nsym = code.rows (); -+ int primpoly = code.primpoly (); -+ int n = args(1).nint_value (); -+ int k = args(2).nint_value (); -+ -+ int m = 1; -+ while (n > (1< __OCTAVE_GALOIS_MAX_M)) -+ { -+ error ("rsdec: invalid values of message and codeword length"); -+ return retval; -+ } -+ -+ if ((n-k) & 1) -+ { -+ error ("rsdec: difference of message and codeword length must be even"); -+ return retval; -+ } -+ -+ int nroots = n-k; -+ galois genpoly; -+ bool have_genpoly = false; -+ bool parity_at_end = true; -+ int fcr = 0; -+ int prim = 0; -+ int iprim; -+ -+ for (int i = 3; i < 6; i++) -+ { -+ if (nargin > i) -+ { -+ if (args(i).is_string ()) -+ { -+ std::string parstr = args(i).string_value (); -+ for (int j = 0; j < (int)parstr.length (); j++) -+ parstr[j] = toupper (parstr[j]); -+ -+ if (!parstr.compare("END")) -+ { -+ parity_at_end = true; -+ } -+ else if (!parstr.compare("BEGINNING")) -+ { -+ parity_at_end = false; -+ } -+ else -+ { -+ error ("rsdec: unrecoginized parrity position"); -+ return retval; -+ } -+ } -+ else -+ { -+ if (args(i).type_id () == octave_galois::static_type_id ()) -+ { -+ if (have_genpoly) -+ { -+ print_usage (); -+ return retval; -+ } -+ genpoly = ((const octave_galois&) args(i).get_rep ()).galois_value (); -+ } -+ else -+ { -+ if (have_genpoly) -+ { -+ if (prim != 0) -+ { -+ print_usage (); -+ return retval; -+ } -+ prim = args(i).nint_value (); -+ } -+ else -+ fcr = args(i).nint_value (); -+ } -+ have_genpoly = true; -+ } -+ } -+ } -+ -+ if (have_genpoly) -+ { -+ if (fcr != 0) -+ { -+ if ((fcr < 1) || (fcr > nn)) -+ { -+ error ("rsdec: invalid first consecutive root of generator polynomial"); -+ return retval; -+ } -+ if ((prim < 1) || (prim > nn)) -+ { -+ error ("rsdec: invalid primitive element of generator polynomial"); -+ return retval; -+ } -+ } -+ else -+ { -+ if (genpoly.cols () > genpoly.rows ()) -+ genpoly = genpoly.transpose (); -+ -+ if (genpoly.cols () != 1) -+ { -+ error ("rsdec: the generator polynomial must be a vector"); -+ return retval; -+ } -+ -+ if (genpoly.primpoly () != primpoly) -+ { -+ error ("rsdec: the generator polynomial must be same galois " -+ "field as the message"); -+ return retval; -+ } -+ -+ if (genpoly.rows () != nroots+1) -+ { -+ error ("rsdec: generator polynomial has incorrect order"); -+ return retval; -+ } -+ -+ // Find the roots of the generator polynomial -+ int count = 0; -+ OCTAVE_LOCAL_BUFFER (int, roots, nroots); -+ for (int j = 0; j <= nn; j++) -+ { -+ // Evaluate generator polynomial at j -+ int val = genpoly(0, 0); -+ int indx = genpoly.index_of (j); -+ for (int i = 0; i 0; i--) -+ code(l, i+nn-n-1) = code(l, i-1); -+ } -+ -+ for (int l = 0; l < nsym; l++) -+ nerr(l) = decode_rs (code, prim, iprim, nroots, fcr, l, parity_at_end); -+ -+ if (nn != n) -+ { -+ if (parity_at_end) -+ for (int l = 0; l < nsym; l++) -+ for (int i = 0; i > n; i--) -+ code(l, i) = code(l, i+nn-n); -+ code.resize (dim_vector (nsym, n), 0); -+ } -+ -+ if (parity_at_end) -+ { -+ for (int l = 0; l < nsym; l++) -+ for (int i = 0; i < k; i++) -+ msg(l, i) = code(l, i); -+ } -+ else -+ { -+ for (int l = 0; l < nsym; l++) -+ for (int i = 0; i < k; i++) -+ msg(l, i) = code(l, nroots+i); -+ } -+ -+ retval(0) = new octave_galois (msg); -+ retval(1) = octave_value (nerr); -+ retval(2) = new octave_galois (code); -+ -+ return retval; -+} -+ -+/* -+%% Test input validation -+%!error rsdec () -+%!error rsdec (1) -+%!error rsdec (1, 2) -+%!error rsdec (1, 2, 3, 4, 5, 6) -+*/ -+ -+// PKG_ADD: autoload ("bchenco", "gf.oct"); -+// PKG_DEL: autoload ("bchenco", "gf.oct", "remove"); -+DEFUN_DLD (bchenco, args, , -+ "-*- texinfo -*-\n\ -+@deftypefn {Loadable Function} {@var{code} =} bchenco (@var{msg}, @var{n}, @var{k})\n\ -+@deftypefnx {Loadable Function} {@var{code} =} bchenco (@var{msg}, @var{n}, @var{k}, @var{g})\n\ -+@deftypefnx {Loadable Function} {@var{code} =} bchenco (@dots{}, @var{parpos})\n\ -+Encodes the message @var{msg} using a [@var{n},@var{k}] BCH coding.\n\ -+The variable @var{msg} is a binary array with @var{k} columns and an\n\ -+arbitrary number of rows. Each row of @var{msg} represents a single symbol\n\ -+to be coded by the BCH coder. The coded message is returned in the binary\n\ -+array @var{code} containing @var{n} columns and the same number of rows as\n\ -+@var{msg}.\n\ -+\n\ -+The use of @code{bchenco} can be seen in the following short example.\n\ -+\n\ -+@example\n\ -+m = 3; n = 2^m -1; k = 4;\n\ -+msg = randint (10,k);\n\ -+code = bchenco (msg, n, k);\n\ -+@end example\n\ -+\n\ -+Valid codes can be found using @code{bchpoly}. In general the codeword\n\ -+length @var{n} should be of the form @code{2^@var{m}-1}, where m is an\n\ -+integer. However, shortened BCH codes can be used such that if\n\ -+@code{[2^@var{m}-1,@var{k}]} is a valid code\n\ -+@code{[2^@var{m}-1-@var{x},@var{k}-@var{x}]}\n is also a valid code using\n\ -+the same generator polynomial.\n\ -+\n\ -+By default the generator polynomial used in the BCH coding is\n\ -+based on the properties of the Galois Field GF(2^@var{m}). This\n\ -+default generator polynomial can be overridden by a polynomial in @var{g}.\n\ -+Suitable generator polynomials can be constructed with @code{bchpoly}.\n\ -+\n\ -+By default the parity symbols are placed at the beginning of the coded\n\ -+message. The variable @var{parpos} controls this positioning and can take\n\ -+the values @code{\"beginning\"} or @code{\"end\"}.\n\ -+@seealso{bchpoly, bchdeco, encode}\n\ -+@end deftypefn") -+{ -+ octave_value retval; -+ int nargin = args.length (); -+ -+ if (nargin < 3 || nargin > 5) -+ { -+ print_usage (); -+ return retval; -+ } -+ -+ Matrix msg = args(0).matrix_value (); -+ int nsym = msg.rows (); -+ int nn = args(1).nint_value (); -+ int k = args(2).nint_value (); -+ -+ int m = 1; -+ while (nn > (1< __OCTAVE_GALOIS_MAX_M)) -+ { -+ error ("bchenco: invalid values of message or codeword length"); -+ return retval; -+ } -+ -+ galois genpoly; -+ bool have_genpoly = false; -+ bool parity_at_end = false; -+ -+ for (int i = 3; i < nargin; i++) -+ { -+ if (args(i).is_string ()) -+ { -+ std::string parstr = args(i).string_value (); -+ for (int j = 0; j < (int)parstr.length (); j++) -+ parstr[j] = toupper (parstr[j]); -+ -+ if (!parstr.compare("END")) -+ { -+ parity_at_end = true; -+ } -+ else if (!parstr.compare("BEGINNING")) -+ { -+ parity_at_end = false; -+ } -+ else -+ { -+ error ("bchenco: unrecoginized parity position"); -+ return retval; -+ } -+ } -+ else -+ { -+ have_genpoly = true; -+ genpoly = galois (args(i).matrix_value (), m); -+ if (genpoly.cols () > genpoly.rows ()) -+ genpoly = genpoly.transpose (); -+ -+ if (genpoly.cols () != 1) -+ { -+ error ("bchenco: the generator polynomial must be a vector"); -+ return retval; -+ } -+ -+ if (genpoly.rows () != nn-k+1) -+ { -+ error ("bchenco: generator polynomial has incorrect order"); -+ return retval; -+ } -+ } -+ } -+ -+ if (!have_genpoly) -+ { -+ // The code below is basically bchpoly.m in C++, so if there is a need -+ // it can be used to rewrite bchpoly as an oct-file... -+ -+ RowVector found (n, 0); -+ found(0) = 1; -+ galois c (1, m, 0, m); -+ c(0, 0) = c.index_of (1); -+ Array cs (dim_vector (1, 1), 1); -+ -+ int nc = 1; -+ -+ // Find the cyclotomic cosets of GF(2^m) -+ while (found.min () == 0) -+ { -+ int idx = n; -+ for (int i = 0; i idx) -+ { -+ c(nc, cs(nc)) = r; -+ found(c.alpha_to (r)-1) = 1; -+ cs(nc) += 1; -+ } -+ nc++; -+ } -+ -+ // Re-use the found vector with 1==not-found !!! -+ found.resize (nc); -+ -+ galois f (1, 0, 0, m); -+ int t = 0; -+ int nf = 0; -+ do -+ { -+ t++; -+ for (int i = 0; i < nc; i++) -+ { -+ if (found(i) == 1) -+ { -+ for (int j = 2*(t-1); j<2*t; j++) -+ { -+ int flag = 0; -+ for (int l = 0; l < cs(i); l++) -+ { -+ if (c(i, l) == j+1) -+ { -+ f.resize (dim_vector (1, nf+cs(i))); -+ for (int ll = 0; ll < cs(i); ll++) -+ f(0, nf+ll) = c(i, ll); -+ found(i) = 0; -+ nf += cs(i); -+ flag = 1; -+ break; -+ } -+ } -+ if (flag) break; -+ } -+ } -+ } -+ } -+ while (nf < nn - k); -+ -+ if (nf != nn - k) -+ { -+ error ("bchenco: can not find valid generator polynomial for parameters"); -+ return retval; -+ } -+ -+ // Create polynomial of right length. -+ genpoly = galois (nf+1, 1, 0, m); -+ -+ genpoly(0, 0) = 1; -+ for (int i = 0; i < nf; i++) -+ { -+ genpoly(i+1, 0) = 1; -+ -+ // Multiply genpoly by @**(root + x) -+ for (int l = i; l > 0; l--) -+ { -+ if (genpoly(l, 0) != 0) -+ genpoly(l, 0) = genpoly(l-1, 0) -+ ^ genpoly.alpha_to (modn (genpoly.index_of (genpoly(l, 0)) + f(0, i), -+ m, n)); -+ else -+ genpoly(l, 0) = genpoly(l-1, 0); -+ } -+ // genpoly(0,0) can never be zero -+ genpoly(0, 0) = genpoly.alpha_to (modn (genpoly.index_of (genpoly(0, 0)) -+ + f(0, i), -+ m, n)); -+ } -+ } -+ -+ // Add space for parity block -+ msg.resize (nsym, nn, 0); -+ -+ // The code below basically finds the parity bits by treating the -+ // message as a polynomial and dividing it by the generator polynomial. -+ // The parity bits are then the remainder of this division. -+ // -+ // This code could just as easily be written as -+ // [ignore par] = gdeconv(gf(msg), gf(genpoly)); -+ // But the code below has the advantage of being 20 times faster :-) -+ -+ if (parity_at_end) -+ { -+ for (int l = 0; l < nsym; l++) -+ { -+ for (int i = 0; i < k; i++) -+ { -+ int feedback = (int)msg(l, i) ^ (int)msg(l, k); -+ if (feedback != 0) -+ { -+ for (int j = 0; j < nn-k-1; j++) -+ if (genpoly(nn-k-j-1, 0) != 0) -+ msg(l, k+j) = (int)msg(l, k+j+1) ^ feedback; -+ else -+ msg(l, k+j) = msg(l, k+j+1); -+ msg(l, nn-1) = genpoly(0, 0) & feedback; -+ } -+ else -+ { -+ for (int j = k; j < nn-1; j++) -+ msg(l, j) = msg(l, j+1); -+ msg(l, nn-1) = 0; -+ } -+ } -+ } -+ } -+ else -+ { -+ for (int l = 0; l < nsym; l++) -+ { -+ for (int i=k; i > 0; i--) -+ msg(l, i+nn-k-1) = msg(l, i-1); -+ for (int i = 0; i= 0; i--) -+ { -+ int feedback = (int)msg(l, nn-k+i) ^ (int)msg(l, nn-k-1); -+ if (feedback != 0) -+ { -+ for (int j = nn - k -1; j > 0; j--) -+ if (genpoly(j, 0) != 0) -+ msg(l, j) = (int)msg(l, j-1) ^ feedback; -+ else -+ msg(l, j) = msg(l, j-1); -+ msg(l, 0) = genpoly(0, 0) & feedback; -+ } -+ else -+ { -+ for (int j = nn - k - 1; j > 0; j--) -+ msg(l, j) = msg(l, j-1); -+ msg(l, 0) = 0; -+ } -+ } -+ } -+ } -+ -+ retval = msg; -+ return retval; -+} -+ -+/* -+%% Test input validation -+%!error bchenco () -+%!error bchenco (1) -+%!error bchenco (1, 2) -+%!error bchenco (1, 2, 3, 4, 5, 6) -+*/ -+ -+// PKG_ADD: autoload ("bchdeco", "gf.oct"); -+// PKG_DEL: autoload ("bchdeco", "gf.oct", "remove"); -+DEFUN_DLD (bchdeco, args, , -+ "-*- texinfo -*-\n\ -+@deftypefn {Loadable Function} {@var{msg} =} bchdeco (@var{code}, @var{k}, @var{t})\n\ -+@deftypefnx {Loadable Function} {@var{msg} =} bchdeco (@var{code}, @var{k}, @var{t}, @var{prim})\n\ -+@deftypefnx {Loadable Function} {@var{msg} =} bchdeco (@dots{}, @var{parpos})\n\ -+@deftypefnx {Loadable Function} {[@var{msg}, @var{err}] =} bchdeco (@dots{})\n\ -+@deftypefnx {Loadable Function} {[@var{msg}, @var{err}, @var{ccode}] =} bchdeco (@dots{})\n\ -+Decodes the coded message @var{code} using a BCH coder. The message length\n\ -+of the coder is defined in variable @var{k}, and the error correction\n\ -+capability of the code is defined in @var{t}.\n\ -+\n\ -+The variable @var{code} is a binary array with @var{n} columns and an\n\ -+arbitrary number of rows. Each row of @var{code} represents a single symbol\n\ -+to be decoded by the BCH coder. The decoded message is returned in the\n\ -+binary array @var{msg} containing @var{k} columns and the same number of\n\ -+rows as @var{code}.\n\ -+\n\ -+The use of @code{bchdeco} can be seen in the following short example.\n\ -+\n\ -+@example\n\ -+m = 3; n = 2^m -1; k = 4; t = 1;\n\ -+msg = randint (10, k);\n\ -+code = bchenco (msg, n, k);\n\ -+noisy = mod (randerr (10,n) + code, 2);\n\ -+[dec, err] = bchdeco (msg, k, t);\n\ -+@end example\n\ -+\n\ -+Valid codes can be found using @code{bchpoly}. In general the codeword\n\ -+length @var{n} should be of the form @code{2^@var{m}-1}, where m is an\n\ -+integer. However, shortened BCH codes can be used such that if\n\ -+@code{[2^@var{m}-1,@var{k}]} is a valid code\n\ -+@code{[2^@var{m}-1-@var{x},@var{k}-@var{x}]}\n is also a valid code using\n\ -+the same generator polynomial.\n\ -+\n\ -+By default the BCH coding is based on the properties of the Galois\n\ -+Field GF(2^@var{m}). The primitive polynomial used in the Galois\n\ -+can be overridden by a primitive polynomial in @var{prim}. Suitable\n\ -+primitive polynomials can be constructed with @code{primpoly}. The form\n\ -+of @var{prim} maybe be either a integer representation of the primitive\n\ -+polynomial as given by @code{primpoly}, or a binary representation that\n\ -+might be constructed like\n\ -+\n\ -+@example\n\ -+m = 3;\n\ -+prim = de2bi (primpoly (m));\n\ -+@end example\n\ -+\n\ -+By default the parity symbols are assumed to be placed at the beginning of\n\ -+the coded message. The variable @var{parpos} controls this positioning and\n\ -+can take the values @code{\"beginning\"} or @code{\"end\"}.\n\ -+@seealso{bchpoly, bchenco, decode, primpoly}\n\ -+@end deftypefn") -+{ -+ octave_value_list retval; -+ int nargin = args.length (); -+ -+ if (nargin < 3 || nargin > 5) -+ { -+ print_usage (); -+ return retval; -+ } -+ -+ Matrix code = args(0).matrix_value (); -+ int nsym = code.rows (); -+ int nn = code.cols (); -+ int k = args(1).nint_value (); -+ int t = args(2).nint_value (); -+ int t2 = t << 1; -+ -+ int m = 1; -+ while (nn > (1< __OCTAVE_GALOIS_MAX_M)) -+ { -+ error ("bchdeco: invalid values of message or codeword length"); -+ return retval; -+ } -+ -+ int prim = 0; // primitve polynomial of zero flags default -+ bool parity_at_end = false; -+ -+ for (int i = 3; i < nargin; i++) -+ { -+ if (args(i).is_string ()) -+ { -+ std::string parstr = args(i).string_value (); -+ for (int j = 0; j < (int)parstr.length (); j++) -+ parstr[j] = toupper (parstr[j]); -+ -+ if (!parstr.compare("END")) -+ { -+ parity_at_end = true; -+ } -+ else if (!parstr.compare("BEGINNING")) -+ { -+ parity_at_end = false; -+ } -+ else -+ { -+ error ("bchdeco: unrecoginized parity position"); -+ return retval; -+ } -+ } -+ else -+ { -+ if (args(i).is_real_scalar ()) -+ prim = args(i).int_value (); -+ else -+ { -+ Matrix tmp = args(i).matrix_value (); -+ -+ if (tmp.cols () > tmp.rows ()) -+ tmp = tmp.transpose (); -+ -+ if (tmp.cols () != 1) -+ { -+ error ("bchdeco: the primitve polynomial must be a scalar " -+ "or a vector"); -+ return retval; -+ } -+ -+ prim = 0; -+ for (int i = 0; i < tmp.rows (); i++) -+ if ((int)tmp(i, 0) & 1) -+ prim |= (1< s (dim_vector(t2+1, 1), 0); -+ bool syn_error = false; -+ -+ for (int i = 1; i <= t2; i++) -+ { -+ for (int j = 0; j < nn; j++) -+ { -+ if (parity_at_end) -+ { -+ if (code(lsym, nn-j-1) != 0) -+ s(i) ^= tables.alpha_to (modn (i*j, m, n)); -+ } -+ else -+ { -+ if (code(lsym, j) != 0) -+ s(i) ^= tables.alpha_to (modn (i*j, m, n)); -+ } -+ } -+ if (s(i) != 0) -+ syn_error = true; /* set error flag if non-zero syndrome */ -+ -+ } -+ -+ if (syn_error) -+ { /* if there are errors, try to correct them */ -+ int q, u; -+ Array d (dim_vector (t2+2, 1)), l(dim_vector (t2+2, 1)), -+ u_lu(dim_vector (t2+2, 1)), reg(dim_vector (t2+2, 1)), -+ elp(dim_vector (t2+2, t2+2)); -+ -+ /* convert syndrome from polynomial form to index form */ -+ for (int i = 1; i <= t2; i++) -+ s(i) = tables.index_of (s(i)); -+ -+ /* -+ * Compute the error location polynomial via the Berlekamp -+ * iterative algorithm. Following the terminology of Lin and -+ * Costello's book : d(u) is the 'mu'th discrepancy, where -+ * u='mu'+1 and 'mu' (the Greek letter!) is the step number -+ * ranging from -1 to 2*t (see L&C), l(u) is the degree of -+ * the elp at that step, and u_l(u) is the difference between -+ * the step number and the degree of the elp. -+ */ -+ /* initialise table entries */ -+ d(0) = 0; /* index form */ -+ d(1) = s(1); /* index form */ -+ elp(0, 0) = 0; /* index form */ -+ elp(1, 0) = 1; /* polynomial form */ -+ for (int i = 1; i < t2; i++) -+ { -+ elp(0, i) = n; /* index form */ -+ elp(1, i) = 0; /* polynomial form */ -+ } -+ l(0) = 0; -+ l(1) = 0; -+ u_lu(0) = -1; -+ u_lu(1) = 0; -+ u = 0; -+ -+ do -+ { -+ u++; -+ if (d(u) == n) -+ { -+ l(u + 1) = l(u); -+ for (int i = 0; i <= l(u); i++) -+ { -+ elp(u + 1, i) = elp(u, i); -+ elp(u, i) = tables.index_of (elp(u, i)); -+ } -+ } -+ else -+ /* -+ * search for words with greatest u_lu(q) for -+ * which d(q)!=0 -+ */ -+ { -+ q = u - 1; -+ while ((d(q) == n) && (q > 0)) -+ q--; -+ /* have found first non-zero d(q) */ -+ if (q > 0) -+ { -+ int j = q; -+ do -+ { -+ j--; -+ if ((d(j) != n) && (u_lu(q) < u_lu(j))) -+ q = j; -+ } -+ while (j > 0); -+ } -+ -+ /* -+ * have now found q such that d(u)!=0 and -+ * u_lu(q) is maximum -+ */ -+ /* store degree of new elp polynomial */ -+ if (l(u) > l(q) + u - q) -+ l(u + 1) = l(u); -+ else -+ l(u + 1) = l(q) + u - q; -+ -+ /* form new elp(x) */ -+ for (int i = 0; i < t2; i++) -+ elp(u + 1, i) = 0; -+ for (int i = 0; i <= l(q); i++) -+ if (elp(q, i) != n) -+ elp(u + 1, i + u - q) = -+ tables.alpha_to (modn ((d(u) + n - d(q) + elp(q, i)), m, n)); -+ for (int i = 0; i <= l(u); i++) -+ { -+ elp(u + 1, i) ^= elp(u, i); -+ elp(u, i) = tables.index_of (elp(u, i)); -+ } -+ } -+ u_lu(u + 1) = u - l(u + 1); -+ -+ /* form (u+1)th discrepancy */ -+ if (u < t2) -+ { -+ /* no discrepancy computed on last iteration */ -+ d(u + 1) = tables.alpha_to (s(u + 1)); -+ -+ for (int i = 1; i <= l(u + 1); i++) -+ if ((s(u + 1 - i) != n) && (elp(u + 1, i) != 0)) -+ d(u + 1) ^= tables.alpha_to (modn (s(u + 1 - i) -+ + tables.index_of (elp(u + 1, i)), -+ m, n)); -+ /* put d(u+1) into index form */ -+ d(u + 1) = tables.index_of (d(u + 1)); -+ } -+ } -+ while ((u < t2) && (l(u + 1) <= t)); -+ -+ u++; -+ if (l(u) <= t) -+ {/* Can correct errors */ -+ int count; -+ Array loc (dim_vector (t+2, 1)); -+ -+ /* put elp into index form */ -+ for (int i = 0; i <= l(u); i++) -+ elp(u, i) = tables.index_of (elp(u, i)); -+ -+ /* Chien search: find roots of the error location polynomial */ -+ for (int i = 1; i <= l(u); i++) -+ reg(i) = elp(u, i); -+ count = 0; -+ for (int i = 1; i <= n; i++) -+ { -+ q = 1; -+ for (int j = 1; j <= l(u); j++) -+ if (reg(j) != n) -+ { -+ reg(j) = modn ((reg(j) + j), m, n); -+ q ^= tables.alpha_to (reg(j)); -+ } -+ if (!q) -+ { /* store root and error -+ * location number indices */ -+ loc(count) = n - i; -+ count++; -+ if (count > l(u)) -+ break; -+ } -+ } -+ -+ if (count == l(u)) -+ { -+ /* no. roots = degree of elp hence <= t errors */ -+ nerr(lsym) = l(u); -+ for (int i = 0; i < l(u); i++) -+ if (parity_at_end) -+ code(lsym, nn-loc(i)-1) = -+ (int)code(lsym, nn-loc(i)-1) ^ 1; -+ else -+ code(lsym, loc(i)) = (int)code(lsym, loc(i)) ^ 1; -+ } -+ else /* elp has degree >t hence cannot solve */ -+ nerr(lsym) = -1; -+ } -+ else -+ nerr(lsym) = -1; -+ } -+ } -+ -+ Matrix msg (nsym, k); -+ if (parity_at_end) -+ { -+ for (int l = 0; l < nsym; l++) -+ for (int i = 0; i < k; i++) -+ msg(l, i) = code(l, i); -+ } -+ else -+ { -+ for (int l = 0; l < nsym; l++) -+ for (int i = 0; i < k; i++) -+ msg(l, i) = code(l, nn-k+i); -+ } -+ -+ retval(0) = octave_value (msg); -+ retval(1) = octave_value (nerr); -+ retval(2) = octave_value (code); -+ return retval; -+} -+ -+/* -+%% Test input validation -+%!error bchdeco () -+%!error bchdeco (1) -+%!error bchdeco (1, 2) -+%!error bchdeco (1, 2, 3, 4, 5, 6) -+*/ -diff -uNr a/src/__gfweight__.cc b/src/__gfweight__.cc ---- a/src/__gfweight__.cc 2015-04-04 12:28:43.938510295 -0400 -+++ b/src/__gfweight__.cc 2018-04-09 13:51:01.973935827 -0400 -@@ -68,7 +68,7 @@ - if (k > 128) - { - octave_stdout << "__gfweight__: this is likely to take a very long time!!\n"; -- flush_octave_stdout (); -+ octave::flush_stdout (); - } - - Array codeword (dim_vector (n, 1), 0); -diff -uNr a/src/__gfweight__.cc~ b/src/__gfweight__.cc~ ---- a/src/__gfweight__.cc~ 1969-12-31 19:00:00.000000000 -0500 -+++ b/src/__gfweight__.cc~ 2015-04-04 12:28:43.938510295 -0400 -@@ -0,0 +1,89 @@ -+//Copyright (C) 2003 David Bateman -+// -+// 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 3 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; if not, see -+// . -+// -+// In addition to the terms of the GPL, you are permitted to link this -+// program with any Open Source program, as defined by the Open Source -+// Initiative (www.opensource.org) -+ -+#include -+ -+static int -+get_weight (const Array& codeword, const Matrix& gen, -+ int weight, int depth, int start, int n, int k) -+{ -+ int retval = weight; -+ -+ for (int i = start; i < k ; i++) -+ { -+ OCTAVE_QUIT; -+ -+ Array new_codeword (codeword); -+ int tmp = 0; -+ for (int j = 0; j < n; j++) -+ if (new_codeword (j) ^= (char)gen(i,j)) -+ tmp++; -+ if (tmp < retval) -+ retval = tmp; -+ if (depth < retval) -+ retval = get_weight (new_codeword, gen, retval, depth+1, i+1, n, k); -+ } -+ return retval; -+} -+ -+DEFUN_DLD (__gfweight__, args, , -+ "-*- texinfo -*-\n\ -+@deftypefn {Loadable Function} {@var{w} =} __gfweight__ (@var{gen})\n\ -+Returns the minimum distance @var{w} of the generator matrix @var{gen}.\n\ -+The codeword length is @var{k}.\n\ -+\n\ -+This is an internal function of @code{gfweight}. You should use\n\ -+@code{gfweight} rather than use this function directly.\n\ -+@seealso{gfweight}\n\ -+@end deftypefn") -+{ -+ -+ if (args.length () != 1) -+ { -+ print_usage (); -+ return octave_value (); -+ } -+ -+ Matrix gen = args(0).matrix_value (); -+ int k = gen.rows (); -+ int n = gen.columns (); -+ -+ if (k > 128) -+ { -+ octave_stdout << "__gfweight__: this is likely to take a very long time!!\n"; -+ flush_octave_stdout (); -+ } -+ -+ Array codeword (dim_vector (n, 1), 0); -+ return octave_value ((double)get_weight (codeword, gen, n - k + 1, 1, -+ 0, n, k)); -+} -+ -+/* -+%% Test input validation -+%!error __gfweight__ () -+%!error __gfweight__ (1, 2) -+*/ -+ -+/* -+;;; Local Variables: *** -+;;; mode: C++ *** -+;;; End: *** -+*/ -diff -uNr a/src/ov-galois.cc b/src/ov-galois.cc ---- a/src/ov-galois.cc 2018-04-09 13:25:42.884981069 -0400 -+++ b/src/ov-galois.cc 2018-04-09 13:27:40.103502409 -0400 -@@ -21,7 +21,7 @@ - #include - - #include --#include -+#include - #include - #include - #include -@@ -328,13 +328,13 @@ - - if (rows () > 0 && columns () > 0) - { -- gripe_implicit_conversion ("Octave:array-as-scalar", -+ warn_implicit_conversion ("Octave:array-as-scalar", - "real matrix", "real scalar"); - - retval = (double) gval (0, 0); - } - else -- gripe_invalid_conversion ("galois", "real scalar"); -+ err_invalid_conversion ("galois", "real scalar"); - - return retval; - } -@@ -348,13 +348,13 @@ - - if (rows () > 0 && columns () > 0) - { -- gripe_implicit_conversion ("Octave:array-as-scalar", -+ warn_implicit_conversion ("Octave:array-as-scalar", - "real matrix", "real scalar"); - - retval = (double) gval (0, 0); - } - else -- gripe_invalid_conversion ("galois", "complex scalar"); -+ err_invalid_conversion ("galois", "complex scalar"); - - return retval; - } -diff -uNr a/src/ov-galois.h b/src/ov-galois.h ---- a/src/ov-galois.h 2018-04-09 13:25:42.872981630 -0400 -+++ b/src/ov-galois.h 2018-04-09 14:09:33.961913785 -0400 -@@ -49,7 +49,6 @@ - #endif - - class octave_value_list; --class tree_walker; - - // Data structures. - -@@ -100,7 +99,7 @@ - - bool is_defined (void) const { return true; } - -- bool is_numeric_type (void) const { return true; } -+ bool isnumeric (void) const { return true; } - - bool is_constant (void) const { return true; } - -@@ -124,7 +123,7 @@ - - bool is_real_matrix (void) const { return false; } - -- bool is_real_type (void) const { return false; } -+ bool isreal (void) const { return false; } - - // FIXME - bool valid_as_scalar_index (void) const { return false; } -diff -uNr a/src/ov-galois.h~ b/src/ov-galois.h~ ---- a/src/ov-galois.h~ 1969-12-31 19:00:00.000000000 -0500 -+++ b/src/ov-galois.h~ 2018-04-09 13:35:21.137945379 -0400 -@@ -0,0 +1,183 @@ -+//Copyright (C) 2003 David Bateman -+// -+// 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 3 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; if not, see -+// . -+// -+// In addition to the terms of the GPL, you are permitted to link this -+// program with any Open Source program, as defined by the Open Source -+// Initiative (www.opensource.org) -+ -+#if !defined (octave_galois_h) -+#define octave_galois_h 1 -+ -+#include -+#include -+ -+#include "galois.h" -+ -+// The keys of the values in the octave map -+#define __GALOIS_PRIMPOLY_STR "prim_poly" -+#define __GALOIS_ORDER_STR "m" -+#define __GALOIS_DATA_STR "x" -+#ifdef GALOIS_DISP_PRIVATES -+#define __GALOIS_LENGTH_STR "n" -+#define __GALOIS_ALPHA_TO_STR "alpha_to" -+#define __GALOIS_INDEX_OF_STR "index_of" -+#endif -+ -+#if !defined (HAVE_OCTAVE_HDF5_ID_TYPE) -+#if defined (HAVE_HDF5) -+typedef hid_t octave_hdf5_id; -+#else -+typedef int octave_hdf5_id; -+#endif -+#endif -+ -+#if ! defined (OV_REP_TYPE) -+# define OV_REP_TYPE octave_base_value -+#endif -+ -+class octave_value_list; -+class tree_walker; -+ -+// Data structures. -+ -+class -+octave_galois : public octave_base_value -+{ -+public: -+ -+ octave_galois (const Matrix& data = Matrix (0, 0), const int _m = 1, -+ const int _primpoly = 0) -+ { gval = galois (data, _m, _primpoly); } -+ -+ octave_galois (const galois& gm) -+ : octave_base_value (), gval (gm) { } -+ octave_galois (const octave_galois& s) -+ : octave_base_value (), gval (s.gval) { } -+ -+ ~octave_galois (void) { }; -+ -+ OV_REP_TYPE *clone (void) const { return new octave_galois (*this); } -+ OV_REP_TYPE *empty_clone (void) const { return new octave_galois (); } -+ -+ octave_value subsref (const std::string &type, -+ const std::list& idx); -+ -+ octave_value_list subsref (const std::string& type, -+ const std::list& idx, int) -+ { return subsref (type, idx); } -+ -+ octave_value do_index_op (const octave_value_list& idx, -+ bool resize_ok); -+ -+ octave_value do_index_op (const octave_value_list& idx) -+ { return do_index_op (idx, 0); } -+ -+ void assign (const octave_value_list& idx, const galois& rhs); -+ -+ dim_vector dims (void) const { return gval.dims (); } -+ -+ octave_value resize (const dim_vector& dv, bool) const; -+ -+ size_t byte_size (void) const { return gval.byte_size (); } -+ -+ octave_value all (int dim = 0) const { return gval.all (dim); } -+ octave_value any (int dim = 0) const { return gval.any (dim); } -+ -+ bool is_matrix_type (void) const { return true; } -+ -+ bool is_defined (void) const { return true; } -+ -+ bool isnumeric (void) const { return true; } -+ -+ bool is_constant (void) const { return true; } -+ -+ bool is_true (void) const; -+ -+ bool is_galois_type (void) const { return true; } -+ -+ bool print_as_scalar (void) const; -+ -+#if defined (HAVE_OCTAVE_BASE_VALUE_PRINT_CONST) -+ void print (std::ostream& os, bool pr_as_read_syntax = false) const; -+#else -+ void print (std::ostream& os, bool pr_as_read_syntax = false); -+#endif -+ -+ void print_raw (std::ostream& os, bool pr_as_read_syntax = false) const; -+ -+ bool print_name_tag (std::ostream& os, const std::string& name) const; -+ -+ void print_info (std::ostream& os, const std::string& prefix) const; -+ -+ bool is_real_matrix (void) const { return false; } -+ -+ bool isreal (void) const { return false; } -+ -+ // FIXME -+ bool valid_as_scalar_index (void) const { return false; } -+ -+ double double_value (bool = false) const; -+ -+ double scalar_value (bool frc_str_conv = false) const -+ { return double_value (frc_str_conv); } -+ -+ Matrix matrix_value (bool = false) const; -+ -+ NDArray array_value (bool = false) const; -+ -+ Complex complex_value (bool = false) const; -+ -+ ComplexMatrix complex_matrix_value (bool = false) const -+ { return ComplexMatrix ( matrix_value ()); } -+ -+ galois galois_value (void) const { return gval; } -+ -+ octave_value_list dotref (const octave_value_list& idx); -+ -+ int m (void) const { return gval.m (); } -+ int primpoly (void) const { return gval.primpoly (); } -+ -+ bool save_ascii (std::ostream& os); -+ -+ bool load_ascii (std::istream& is); -+ -+ bool save_binary (std::ostream& os, bool& save_as_floats); -+ -+ bool load_binary (std::istream& is, bool swap, -+ oct_mach_info::float_format fmt); -+ -+ bool save_hdf5 (octave_hdf5_id loc_id, const char *name, bool save_as_floats); -+ -+ bool load_hdf5 (octave_hdf5_id loc_id, const char *name); -+ -+private: -+ // The array used to managed the Galios Field data -+ galois gval; -+ -+#if defined (DECLARE_OCTAVE_ALLOCATOR) -+ DECLARE_OCTAVE_ALLOCATOR -+#endif -+ -+ DECLARE_OV_TYPEID_FUNCTIONS_AND_DATA -+}; -+ -+#endif -+ -+/* -+;;; Local Variables: *** -+;;; mode: C++ *** -+;;; End: *** -+*/ diff -r 74fe6c154dee -r aad0cc165a27 src/of-communications-6-deprecated.patch --- a/src/of-communications-6-deprecated.patch Wed Jan 01 18:52:37 2020 -0500 +++ /dev/null Thu Jan 01 00:00:00 1970 +0000 @@ -1,438 +0,0 @@ -diff -r 91279e205a70 src/galois.h ---- a/src/galois.h Sun Jan 20 19:13:40 2019 +0100 -+++ b/src/galois.h Sun Jan 20 19:55:05 2019 +0100 -@@ -22,6 +22,7 @@ - #define octave_galois_int_h 1 - - #include -+#include - - #include "galoisfield.h" - #include "base-lu.h" -@@ -150,12 +151,12 @@ - pivot_type ptype; - }; - --void install_gm_gm_ops (void); --void install_gm_m_ops (void); --void install_m_gm_ops (void); --void install_gm_s_ops (void); --void install_s_gm_ops (void); --void install_fil_gm_ops (void); -+void install_gm_gm_ops (octave::type_info& ti); -+void install_gm_m_ops (octave::type_info& ti); -+void install_m_gm_ops (octave::type_info& ti); -+void install_gm_s_ops (octave::type_info& ti); -+void install_s_gm_ops (octave::type_info& ti); -+void install_fil_gm_ops (octave::type_info& ti); - - galois elem_pow (const galois& a, const galois& b); - galois elem_pow (const galois& a, const Matrix& b); -diff -r 91279e205a70 src/gf.cc ---- a/src/gf.cc Sun Jan 20 19:13:40 2019 +0100 -+++ b/src/gf.cc Sun Jan 20 19:55:05 2019 +0100 -@@ -116,13 +116,16 @@ - if (!galois_type_loaded) - { - octave_galois::register_type (); -- install_gm_gm_ops (); -- install_m_gm_ops (); -- install_gm_m_ops (); -- install_s_gm_ops (); -- install_gm_s_ops (); -+ interp.mlock (); -+ -+ octave::type_info& ti = interp.get_type_info (); -+ -+ install_gm_gm_ops (ti); -+ install_m_gm_ops (ti); -+ install_gm_m_ops (ti); -+ install_s_gm_ops (ti); -+ install_gm_s_ops (ti); - galois_type_loaded = true; -- interp.mlock (); - } - - retval = new octave_galois (data, m, primpoly); -diff -r 91279e205a70 src/op-gm-gm.cc ---- a/src/op-gm-gm.cc Sun Jan 20 19:13:40 2019 +0100 -+++ b/src/op-gm-gm.cc Sun Jan 20 19:55:05 2019 +0100 -@@ -18,7 +18,7 @@ - // program with any Open Source program, as defined by the Open Source - // Initiative (www.opensource.org) - --#include -+#include - #include - - #include "galois.h" -@@ -86,36 +86,36 @@ - DEFASSIGNOP_FN (assign, galois, galois, assign) - - void --install_gm_gm_ops (void) -+install_gm_gm_ops (octave::type_info& ti) - { -- INSTALL_UNOP (op_not, octave_galois, not); -- INSTALL_UNOP (op_uminus, octave_galois, uminus); -- INSTALL_UNOP (op_uplus, octave_galois, uplus); -- INSTALL_UNOP (op_transpose, octave_galois, transpose); -- INSTALL_UNOP (op_hermitian, octave_galois, transpose); -+ INSTALL_UNOP_TI (ti, op_not, octave_galois, not); -+ INSTALL_UNOP_TI (ti, op_uminus, octave_galois, uminus); -+ INSTALL_UNOP_TI (ti, op_uplus, octave_galois, uplus); -+ INSTALL_UNOP_TI (ti, op_transpose, octave_galois, transpose); -+ INSTALL_UNOP_TI (ti, op_hermitian, octave_galois, transpose); - -- INSTALL_BINOP (op_add, octave_galois, octave_galois, add); -- INSTALL_BINOP (op_sub, octave_galois, octave_galois, sub); -- INSTALL_BINOP (op_mul, octave_galois, octave_galois, mul); -- INSTALL_BINOP (op_div, octave_galois, octave_galois, div); -- INSTALL_BINOP (op_pow, octave_galois, octave_galois, pow); -- INSTALL_BINOP (op_ldiv, octave_galois, octave_galois, ldiv); -- INSTALL_BINOP (op_lt, octave_galois, octave_galois, lt); -- INSTALL_BINOP (op_le, octave_galois, octave_galois, le); -- INSTALL_BINOP (op_eq, octave_galois, octave_galois, eq); -- INSTALL_BINOP (op_ge, octave_galois, octave_galois, ge); -- INSTALL_BINOP (op_gt, octave_galois, octave_galois, gt); -- INSTALL_BINOP (op_ne, octave_galois, octave_galois, ne); -- INSTALL_BINOP (op_el_mul, octave_galois, octave_galois, el_mul); -- INSTALL_BINOP (op_el_div, octave_galois, octave_galois, el_div); -- INSTALL_BINOP (op_el_pow, octave_galois, octave_galois, el_pow); -- INSTALL_BINOP (op_el_ldiv, octave_galois, octave_galois, el_ldiv); -- INSTALL_BINOP (op_el_and, octave_galois, octave_galois, el_and); -- INSTALL_BINOP (op_el_or, octave_galois, octave_galois, el_or); -+ INSTALL_BINOP_TI (ti, op_add, octave_galois, octave_galois, add); -+ INSTALL_BINOP_TI (ti, op_sub, octave_galois, octave_galois, sub); -+ INSTALL_BINOP_TI (ti, op_mul, octave_galois, octave_galois, mul); -+ INSTALL_BINOP_TI (ti, op_div, octave_galois, octave_galois, div); -+ INSTALL_BINOP_TI (ti, op_pow, octave_galois, octave_galois, pow); -+ INSTALL_BINOP_TI (ti, op_ldiv, octave_galois, octave_galois, ldiv); -+ INSTALL_BINOP_TI (ti, op_lt, octave_galois, octave_galois, lt); -+ INSTALL_BINOP_TI (ti, op_le, octave_galois, octave_galois, le); -+ INSTALL_BINOP_TI (ti, op_eq, octave_galois, octave_galois, eq); -+ INSTALL_BINOP_TI (ti, op_ge, octave_galois, octave_galois, ge); -+ INSTALL_BINOP_TI (ti, op_gt, octave_galois, octave_galois, gt); -+ INSTALL_BINOP_TI (ti, op_ne, octave_galois, octave_galois, ne); -+ INSTALL_BINOP_TI (ti, op_el_mul, octave_galois, octave_galois, el_mul); -+ INSTALL_BINOP_TI (ti, op_el_div, octave_galois, octave_galois, el_div); -+ INSTALL_BINOP_TI (ti, op_el_pow, octave_galois, octave_galois, el_pow); -+ INSTALL_BINOP_TI (ti, op_el_ldiv, octave_galois, octave_galois, el_ldiv); -+ INSTALL_BINOP_TI (ti, op_el_and, octave_galois, octave_galois, el_and); -+ INSTALL_BINOP_TI (ti, op_el_or, octave_galois, octave_galois, el_or); - -- INSTALL_G_CATOP (octave_galois, octave_galois, gm_gm); -+ INSTALL_CATOP_TI (ti, octave_galois, octave_galois, gm_gm); - -- INSTALL_ASSIGNOP (op_asn_eq, octave_galois, octave_galois, assign); -+ INSTALL_ASSIGNOP_TI (ti, op_asn_eq, octave_galois, octave_galois, assign); - } - - /* -diff -r 91279e205a70 src/op-gm-m.cc ---- a/src/op-gm-m.cc Sun Jan 20 19:13:40 2019 +0100 -+++ b/src/op-gm-m.cc Sun Jan 20 19:55:05 2019 +0100 -@@ -75,31 +75,31 @@ - } - - void --install_gm_m_ops (void) -+install_gm_m_ops (octave::type_info& ti) - { -- INSTALL_BINOP (op_add, octave_galois, octave_matrix, add); -- INSTALL_BINOP (op_sub, octave_galois, octave_matrix, sub); -- INSTALL_BINOP (op_mul, octave_galois, octave_matrix, mul); -- INSTALL_BINOP (op_div, octave_galois, octave_matrix, div); -- INSTALL_BINOP (op_pow, octave_galois, octave_matrix, pow); -- INSTALL_BINOP (op_ldiv, octave_galois, octave_matrix, ldiv); -- INSTALL_BINOP (op_lt, octave_galois, octave_matrix, lt); -- INSTALL_BINOP (op_le, octave_galois, octave_matrix, le); -- INSTALL_BINOP (op_eq, octave_galois, octave_matrix, eq); -- INSTALL_BINOP (op_ge, octave_galois, octave_matrix, ge); -- INSTALL_BINOP (op_gt, octave_galois, octave_matrix, gt); -- INSTALL_BINOP (op_ne, octave_galois, octave_matrix, ne); -- INSTALL_BINOP (op_el_mul, octave_galois, octave_matrix, el_mul); -- INSTALL_BINOP (op_el_div, octave_galois, octave_matrix, el_div); -- INSTALL_BINOP (op_el_pow, octave_galois, octave_matrix, el_pow); -- INSTALL_BINOP (op_el_ldiv, octave_galois, octave_matrix, el_ldiv); -- INSTALL_BINOP (op_el_and, octave_galois, octave_matrix, el_and); -- INSTALL_BINOP (op_el_or, octave_galois, octave_matrix, el_or); -+ INSTALL_BINOP_TI (ti, op_add, octave_galois, octave_matrix, add); -+ INSTALL_BINOP_TI (ti, op_sub, octave_galois, octave_matrix, sub); -+ INSTALL_BINOP_TI (ti, op_mul, octave_galois, octave_matrix, mul); -+ INSTALL_BINOP_TI (ti, op_div, octave_galois, octave_matrix, div); -+ INSTALL_BINOP_TI (ti, op_pow, octave_galois, octave_matrix, pow); -+ INSTALL_BINOP_TI (ti, op_ldiv, octave_galois, octave_matrix, ldiv); -+ INSTALL_BINOP_TI (ti, op_lt, octave_galois, octave_matrix, lt); -+ INSTALL_BINOP_TI (ti, op_le, octave_galois, octave_matrix, le); -+ INSTALL_BINOP_TI (ti, op_eq, octave_galois, octave_matrix, eq); -+ INSTALL_BINOP_TI (ti, op_ge, octave_galois, octave_matrix, ge); -+ INSTALL_BINOP_TI (ti, op_gt, octave_galois, octave_matrix, gt); -+ INSTALL_BINOP_TI (ti, op_ne, octave_galois, octave_matrix, ne); -+ INSTALL_BINOP_TI (ti, op_el_mul, octave_galois, octave_matrix, el_mul); -+ INSTALL_BINOP_TI (ti, op_el_div, octave_galois, octave_matrix, el_div); -+ INSTALL_BINOP_TI (ti, op_el_pow, octave_galois, octave_matrix, el_pow); -+ INSTALL_BINOP_TI (ti, op_el_ldiv, octave_galois, octave_matrix, el_ldiv); -+ INSTALL_BINOP_TI (ti, op_el_and, octave_galois, octave_matrix, el_and); -+ INSTALL_BINOP_TI (ti, op_el_or, octave_galois, octave_matrix, el_or); - -- INSTALL_G_CATOP (octave_galois, octave_matrix, gm_m); -+ INSTALL_CATOP_TI (ti, octave_galois, octave_matrix, gm_m); - -- INSTALL_ASSIGNOP (op_asn_eq, octave_galois, octave_matrix, assign); -- INSTALL_ASSIGNCONV (octave_base_value, octave_galois, octave_matrix); -+ INSTALL_ASSIGNOP_TI (ti, op_asn_eq, octave_galois, octave_matrix, assign); -+ INSTALL_ASSIGNCONV_TI (ti, octave_base_value, octave_galois, octave_matrix); - } - - /* -diff -r 91279e205a70 src/op-gm-s.cc ---- a/src/op-gm-s.cc Sun Jan 20 19:13:40 2019 +0100 -+++ b/src/op-gm-s.cc Sun Jan 20 19:55:05 2019 +0100 -@@ -79,30 +79,30 @@ - } - - void --install_gm_s_ops (void) -+install_gm_s_ops (octave::type_info& ti) - { -- INSTALL_BINOP (op_add, octave_galois, octave_scalar, add); -- INSTALL_BINOP (op_sub, octave_galois, octave_scalar, sub); -- INSTALL_BINOP (op_mul, octave_galois, octave_scalar, mul); -- INSTALL_BINOP (op_div, octave_galois, octave_scalar, div); -- INSTALL_BINOP (op_pow, octave_galois, octave_scalar, pow); -- INSTALL_BINOP (op_ldiv, octave_galois, octave_scalar, ldiv); -- INSTALL_BINOP (op_lt, octave_galois, octave_scalar, lt); -- INSTALL_BINOP (op_le, octave_galois, octave_scalar, le); -- INSTALL_BINOP (op_eq, octave_galois, octave_scalar, eq); -- INSTALL_BINOP (op_ge, octave_galois, octave_scalar, ge); -- INSTALL_BINOP (op_gt, octave_galois, octave_scalar, gt); -- INSTALL_BINOP (op_ne, octave_galois, octave_scalar, ne); -- INSTALL_BINOP (op_el_mul, octave_galois, octave_scalar, el_mul); -- INSTALL_BINOP (op_el_div, octave_galois, octave_scalar, el_div); -- INSTALL_BINOP (op_el_pow, octave_galois, octave_scalar, el_pow); -- INSTALL_BINOP (op_el_ldiv, octave_galois, octave_scalar, el_ldiv); -- INSTALL_BINOP (op_el_and, octave_galois, octave_scalar, el_and); -- INSTALL_BINOP (op_el_or, octave_galois, octave_scalar, el_or); -+ INSTALL_BINOP_TI (ti, op_add, octave_galois, octave_scalar, add); -+ INSTALL_BINOP_TI (ti, op_sub, octave_galois, octave_scalar, sub); -+ INSTALL_BINOP_TI (ti, op_mul, octave_galois, octave_scalar, mul); -+ INSTALL_BINOP_TI (ti, op_div, octave_galois, octave_scalar, div); -+ INSTALL_BINOP_TI (ti, op_pow, octave_galois, octave_scalar, pow); -+ INSTALL_BINOP_TI (ti, op_ldiv, octave_galois, octave_scalar, ldiv); -+ INSTALL_BINOP_TI (ti, op_lt, octave_galois, octave_scalar, lt); -+ INSTALL_BINOP_TI (ti, op_le, octave_galois, octave_scalar, le); -+ INSTALL_BINOP_TI (ti, op_eq, octave_galois, octave_scalar, eq); -+ INSTALL_BINOP_TI (ti, op_ge, octave_galois, octave_scalar, ge); -+ INSTALL_BINOP_TI (ti, op_gt, octave_galois, octave_scalar, gt); -+ INSTALL_BINOP_TI (ti, op_ne, octave_galois, octave_scalar, ne); -+ INSTALL_BINOP_TI (ti, op_el_mul, octave_galois, octave_scalar, el_mul); -+ INSTALL_BINOP_TI (ti, op_el_div, octave_galois, octave_scalar, el_div); -+ INSTALL_BINOP_TI (ti, op_el_pow, octave_galois, octave_scalar, el_pow); -+ INSTALL_BINOP_TI (ti, op_el_ldiv, octave_galois, octave_scalar, el_ldiv); -+ INSTALL_BINOP_TI (ti, op_el_and, octave_galois, octave_scalar, el_and); -+ INSTALL_BINOP_TI (ti, op_el_or, octave_galois, octave_scalar, el_or); - -- INSTALL_G_CATOP (octave_galois, octave_scalar, gm_s); -+ INSTALL_CATOP_TI (ti, octave_galois, octave_scalar, gm_s); - -- INSTALL_ASSIGNOP (op_asn_eq, octave_galois, octave_scalar, assign); -+ INSTALL_ASSIGNOP_TI (ti, op_asn_eq, octave_galois, octave_scalar, assign); - } - - /* -diff -r 91279e205a70 src/op-m-gm.cc ---- a/src/op-m-gm.cc Sun Jan 20 19:13:40 2019 +0100 -+++ b/src/op-m-gm.cc Sun Jan 20 19:55:05 2019 +0100 -@@ -76,30 +76,30 @@ - DEFASSIGNOP_FN (assign, matrix, galois, assign) - - void --install_m_gm_ops (void) -+install_m_gm_ops (octave::type_info& ti) - { -- INSTALL_BINOP (op_add, octave_matrix, octave_galois, add); -- INSTALL_BINOP (op_sub, octave_matrix, octave_galois, sub); -- INSTALL_BINOP (op_mul, octave_matrix, octave_galois, mul); -- INSTALL_BINOP (op_div, octave_matrix, octave_galois, div); -- INSTALL_BINOP (op_pow, octave_matrix, octave_galois, pow); -- INSTALL_BINOP (op_ldiv, octave_matrix, octave_galois, ldiv); -- INSTALL_BINOP (op_lt, octave_matrix, octave_galois, lt); -- INSTALL_BINOP (op_le, octave_matrix, octave_galois, le); -- INSTALL_BINOP (op_eq, octave_matrix, octave_galois, eq); -- INSTALL_BINOP (op_ge, octave_matrix, octave_galois, ge); -- INSTALL_BINOP (op_gt, octave_matrix, octave_galois, gt); -- INSTALL_BINOP (op_ne, octave_matrix, octave_galois, ne); -- INSTALL_BINOP (op_el_mul, octave_matrix, octave_galois, el_mul); -- INSTALL_BINOP (op_el_div, octave_matrix, octave_galois, el_div); -- INSTALL_BINOP (op_el_pow, octave_matrix, octave_galois, el_pow); -- INSTALL_BINOP (op_el_ldiv, octave_matrix, octave_galois, el_ldiv); -- INSTALL_BINOP (op_el_and, octave_matrix, octave_galois, el_and); -- INSTALL_BINOP (op_el_or, octave_matrix, octave_galois, el_or); -+ INSTALL_BINOP_TI (ti, op_add, octave_matrix, octave_galois, add); -+ INSTALL_BINOP_TI (ti, op_sub, octave_matrix, octave_galois, sub); -+ INSTALL_BINOP_TI (ti, op_mul, octave_matrix, octave_galois, mul); -+ INSTALL_BINOP_TI (ti, op_div, octave_matrix, octave_galois, div); -+ INSTALL_BINOP_TI (ti, op_pow, octave_matrix, octave_galois, pow); -+ INSTALL_BINOP_TI (ti, op_ldiv, octave_matrix, octave_galois, ldiv); -+ INSTALL_BINOP_TI (ti, op_lt, octave_matrix, octave_galois, lt); -+ INSTALL_BINOP_TI (ti, op_le, octave_matrix, octave_galois, le); -+ INSTALL_BINOP_TI (ti, op_eq, octave_matrix, octave_galois, eq); -+ INSTALL_BINOP_TI (ti, op_ge, octave_matrix, octave_galois, ge); -+ INSTALL_BINOP_TI (ti, op_gt, octave_matrix, octave_galois, gt); -+ INSTALL_BINOP_TI (ti, op_ne, octave_matrix, octave_galois, ne); -+ INSTALL_BINOP_TI (ti, op_el_mul, octave_matrix, octave_galois, el_mul); -+ INSTALL_BINOP_TI (ti, op_el_div, octave_matrix, octave_galois, el_div); -+ INSTALL_BINOP_TI (ti, op_el_pow, octave_matrix, octave_galois, el_pow); -+ INSTALL_BINOP_TI (ti, op_el_ldiv, octave_matrix, octave_galois, el_ldiv); -+ INSTALL_BINOP_TI (ti, op_el_and, octave_matrix, octave_galois, el_and); -+ INSTALL_BINOP_TI (ti, op_el_or, octave_matrix, octave_galois, el_or); - -- INSTALL_G_CATOP (octave_matrix, octave_galois, m_gm); -+ INSTALL_CATOP_TI (ti, octave_matrix, octave_galois, m_gm); - -- INSTALL_ASSIGNOP (op_asn_eq, octave_matrix, octave_galois, assign); -+ INSTALL_ASSIGNOP_TI (ti, op_asn_eq, octave_matrix, octave_galois, assign); - //INSTALL_ASSIGNCONV (octave_base_value, octave_matrix, octave_galois); - } - -diff -r 91279e205a70 src/op-s-gm.cc ---- a/src/op-s-gm.cc Sun Jan 20 19:13:40 2019 +0100 -+++ b/src/op-s-gm.cc Sun Jan 20 19:55:05 2019 +0100 -@@ -83,30 +83,30 @@ - } - - void --install_s_gm_ops (void) -+install_s_gm_ops (octave::type_info& ti) - { -- INSTALL_BINOP (op_add, octave_scalar, octave_galois, add); -- INSTALL_BINOP (op_sub, octave_scalar, octave_galois, sub); -- INSTALL_BINOP (op_mul, octave_scalar, octave_galois, mul); -- INSTALL_BINOP (op_div, octave_scalar, octave_galois, div); -- INSTALL_BINOP (op_pow, octave_scalar, octave_galois, pow); -- INSTALL_BINOP (op_ldiv, octave_scalar, octave_galois, ldiv); -- INSTALL_BINOP (op_lt, octave_scalar, octave_galois, lt); -- INSTALL_BINOP (op_le, octave_scalar, octave_galois, le); -- INSTALL_BINOP (op_eq, octave_scalar, octave_galois, eq); -- INSTALL_BINOP (op_ge, octave_scalar, octave_galois, ge); -- INSTALL_BINOP (op_gt, octave_scalar, octave_galois, gt); -- INSTALL_BINOP (op_ne, octave_scalar, octave_galois, ne); -- INSTALL_BINOP (op_el_mul, octave_scalar, octave_galois, el_mul); -- INSTALL_BINOP (op_el_div, octave_scalar, octave_galois, el_div); -- INSTALL_BINOP (op_el_pow, octave_scalar, octave_galois, el_pow); -- INSTALL_BINOP (op_el_ldiv, octave_scalar, octave_galois, el_ldiv); -- INSTALL_BINOP (op_el_and, octave_scalar, octave_galois, el_and); -- INSTALL_BINOP (op_el_or, octave_scalar, octave_galois, el_or); -+ INSTALL_BINOP_TI (ti, op_add, octave_scalar, octave_galois, add); -+ INSTALL_BINOP_TI (ti, op_sub, octave_scalar, octave_galois, sub); -+ INSTALL_BINOP_TI (ti, op_mul, octave_scalar, octave_galois, mul); -+ INSTALL_BINOP_TI (ti, op_div, octave_scalar, octave_galois, div); -+ INSTALL_BINOP_TI (ti, op_pow, octave_scalar, octave_galois, pow); -+ INSTALL_BINOP_TI (ti, op_ldiv, octave_scalar, octave_galois, ldiv); -+ INSTALL_BINOP_TI (ti, op_lt, octave_scalar, octave_galois, lt); -+ INSTALL_BINOP_TI (ti, op_le, octave_scalar, octave_galois, le); -+ INSTALL_BINOP_TI (ti, op_eq, octave_scalar, octave_galois, eq); -+ INSTALL_BINOP_TI (ti, op_ge, octave_scalar, octave_galois, ge); -+ INSTALL_BINOP_TI (ti, op_gt, octave_scalar, octave_galois, gt); -+ INSTALL_BINOP_TI (ti, op_ne, octave_scalar, octave_galois, ne); -+ INSTALL_BINOP_TI (ti, op_el_mul, octave_scalar, octave_galois, el_mul); -+ INSTALL_BINOP_TI (ti, op_el_div, octave_scalar, octave_galois, el_div); -+ INSTALL_BINOP_TI (ti, op_el_pow, octave_scalar, octave_galois, el_pow); -+ INSTALL_BINOP_TI (ti, op_el_ldiv, octave_scalar, octave_galois, el_ldiv); -+ INSTALL_BINOP_TI (ti, op_el_and, octave_scalar, octave_galois, el_and); -+ INSTALL_BINOP_TI (ti, op_el_or, octave_scalar, octave_galois, el_or); - -- INSTALL_G_CATOP (octave_scalar, octave_galois, s_gm); -+ INSTALL_CATOP_TI (ti, octave_scalar, octave_galois, s_gm); - -- INSTALL_ASSIGNCONV (octave_scalar, octave_galois, octave_galois); -+ INSTALL_ASSIGNCONV_TI (ti, octave_scalar, octave_galois, octave_galois); - } - - /* -diff -r 91279e205a70 src/ov-galois.cc ---- a/src/ov-galois.cc Sun Jan 20 19:13:40 2019 +0100 -+++ b/src/ov-galois.cc Sun Jan 20 19:55:05 2019 +0100 -@@ -24,8 +24,8 @@ - #include - #include - #include --#include - #include -+#include - #include - - -@@ -535,7 +535,7 @@ - - bool - octave_galois::load_binary (std::istream& is, bool swap, -- oct_mach_info::float_format fmt) -+ octave::mach_info::float_format fmt) - { - char mord; - int32_t prim, mdims; -diff -r 91279e205a70 src/ov-galois.h ---- a/src/ov-galois.h Sun Jan 20 19:13:40 2019 +0100 -+++ b/src/ov-galois.h Sun Jan 20 19:55:05 2019 +0100 -@@ -156,7 +156,7 @@ - bool save_binary (std::ostream& os, bool& save_as_floats); - - bool load_binary (std::istream& is, bool swap, -- oct_mach_info::float_format fmt); -+ octave::mach_info::float_format fmt); - - bool save_hdf5 (octave_hdf5_id loc_id, const char *name, bool save_as_floats); - -diff -r 91279e205a70 src/primpoly.cc ---- a/src/primpoly.cc Sun Jan 20 19:13:40 2019 +0100 -+++ b/src/primpoly.cc Sun Jan 20 19:55:05 2019 +0100 -@@ -216,8 +216,8 @@ - for (int i = (1<