changeset 5256:aad0cc165a27

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
author John Donoghue <john.donoghue@ieee.org>
date Mon, 06 Jan 2020 06:19:36 -0500
parents 74fe6c154dee
children 1927e90b6f63
files dist-files.mk 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 src/of-communications.mk
diffstat 8 files changed, 2 insertions(+), 7176 deletions(-) [+]
line wrap: on
line diff
--- 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 \
--- 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 <octave/config.h>
-+#incude <octave/ov.h>
-+])
-+AC_CHECK_HEADER([octave/ls-oct-text.h], [PKG_CPPFLAGS="$PKG_CPPFLAGS -DHAVE_OCTAVE_TEXT_H=1"], [], 
-+[
-+#include <octave/config.h>
-+#include <octave/ov.h>
-+])
-+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 <octave/ov.h>
- #include <octave/pr-output.h>
- 
--#include <octave/ls-oct-ascii.h>
-+
-+#if defined(HAVE_OCTAVE_TEXT_H)
-+#  include <octave/ls-oct-text.h>
-+#else
-+#  include <octave/ls-oct-ascii.h>
-+#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;
- 
--- 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
-+<http://www.gnu.org/licenses/>.
-+
-+*/
-+
-+#ifdef HAVE_CONFIG_H
-+#include <config.h>
-+#endif
-+
-+#include "base-lu.h"
-+
-+template <class lu_type>
-+base_lu<lu_type>::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 <class lu_type>
-+bool
-+base_lu <lu_type> :: packed (void) const
-+{
-+  return l_fact.dims () == dim_vector ();
-+}
-+
-+template <class lu_type>
-+void
-+base_lu <lu_type> :: unpack (void)
-+{
-+  if (packed ())
-+    {
-+      l_fact = L ();
-+      a_fact = U (); // FIXME: sub-optimal
-+      ipvt = getp ();
-+    }
-+}
-+
-+template <class lu_type>
-+lu_type
-+base_lu <lu_type> :: 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 <class lu_type>
-+lu_type
-+base_lu <lu_type> :: 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 <class lu_type>
-+lu_type
-+base_lu <lu_type> :: Y (void) const
-+{
-+  if (! packed ())
-+    (*current_liboctave_error_handler)
-+      ("lu: Y () not implemented for unpacked form");
-+  return a_fact;
-+}
-+
-+template <class lu_type>
-+Array<octave_idx_type>
-+base_lu <lu_type> :: getp (void) const
-+{
-+  if (packed ())
-+    {
-+      octave_idx_type a_nr = a_fact.rows ();
-+
-+      Array<octave_idx_type> 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 <class lu_type>
-+PermMatrix
-+base_lu <lu_type> :: P (void) const
-+{
-+  return PermMatrix (getp (), false);
-+}
-+
-+template <class lu_type>
-+ColumnVector
-+base_lu <lu_type> :: P_vec (void) const
-+{
-+  octave_idx_type a_nr = a_fact.rows ();
-+
-+  ColumnVector p (a_nr);
-+
-+  Array<octave_idx_type> pvt = getp ();
-+
-+  for (octave_idx_type i = 0; i < a_nr; i++)
-+    p.xelem (i) = static_cast<double> (pvt.xelem (i) + 1);
-+
-+  return p;
-+}
-+
-+template <class lu_type>
-+bool
-+base_lu<lu_type>::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
-+<http://www.gnu.org/licenses/>.
-+
-+*/
-+
-+#if !defined (octave_base_lu_h)
-+#define octave_base_lu_h 1
-+
-+#include "MArray.h"
-+#include "dColVector.h"
-+#include "PermMatrix.h"
-+
-+template <class lu_type>
-+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<octave_idx_type> getp (void) const;
-+
-+  lu_type a_fact;
-+  lu_type l_fact;
-+
-+  Array<octave_idx_type> 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<t1> (a1);                \
-+  t2 v2 = dynamic_cast<t2> (a2)
-+#endif
-+
-+#if ! defined (CAST_UNOP_ARG)
-+#  define CAST_UNOP_ARG(t)                      \
-+  t v = dynamic_cast<t> (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<octave_idx_type>& 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 <octave/base-lu.cc>
-+#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 <octave/config.h>
--#include <octave/base-lu.h>
- #include <octave/mx-base.h>
- 
- #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 <octave/byte-swap.h>
- #include <octave/gripes.h>
- #include <octave/lo-ieee.h>
--#include <octave/oct-hdf5.h>
- #include <octave/oct-locbuf.h>
- #include <octave/oct-obj.h>
- #include <octave/ov.h>
-@@ -35,6 +34,22 @@
- #include "galois.h"
- #include "ov-galois.h"
- 
-+#if defined (HAVE_HDF5)
-+#  if defined (HAVE_HDF5_H)
-+#    include <hdf5.h>
-+#  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 <octave/oct.h>
- 
--#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<int>
- get_errs (const int& nmin, const int& nmax, const int &nerrs)
--- 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 <octave/config.h>
-        #include <octave/ov-base.h>
- 
- int
-@@ -2502,7 +2501,6 @@
-    cat confdefs.h - <<_ACEOF >conftest.$ac_ext
- /* end confdefs.h.  */
- 
--       #include <octave/config.h>
-        #include <octave/ov-base.h>
-        class foo : public octave_base_value
-        {
-@@ -2561,7 +2559,6 @@
- /* end confdefs.h.  */
- 
-        #include <iostream>
--       #include <octave/config.h>
-        #include <octave/ov-base.h>
- 
- 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 <octave/config.h>
-        #include <octave/ov-base.h>
-        ]], [[
-        octave_hdf5_id x;
-@@ -56,7 +55,6 @@
-    CXXFLAGS="$CXXFLAGS $comm_CXXFLAGS"
-    AC_LANG_PUSH(C++)
-    AC_COMPILE_IFELSE([AC_LANG_PROGRAM([[
--       #include <octave/config.h>
-        #include <octave/ov-base.h>
-        class foo : public octave_base_value
-        {
-@@ -85,7 +83,6 @@
-    AC_LANG_PUSH(C++)
-    AC_COMPILE_IFELSE([AC_LANG_PROGRAM([[
-        #include <iostream>
--       #include <octave/config.h>
-        #include <octave/ov-base.h>
-        ]], [[
-        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 <octave/config.h>
- #include <octave/error.h>
- #include <octave/gripes.h>
- #include <octave/mx-op-defs.h>
-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 <octave/config.h>
- #include <octave/mx-base.h>
- 
- #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 <octave/config.h>
- #include <octave/defun-dld.h>
- #include <octave/gripes.h>
- #include <octave/oct-locbuf.h>
-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 <octave/config.h>
- #include <octave/oct-obj.h>
- #include <octave/ops.h>
- 
-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 <octave/config.h>
- #include <octave/ops.h>
- #include <octave/ov-re-mat.h>
- 
-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 <octave/config.h>
- #include <octave/ops.h>
- #include <octave/ov-scalar.h>
- 
-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 <octave/config.h>
- #include <octave/ops.h>
- #include <octave/ov-re-mat.h>
- 
-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 <octave/config.h>
- #include <octave/ops.h>
- #include <octave/ov-scalar.h>
- 
-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 <iostream>
- 
--#include <octave/config.h>
- #include <octave/byte-swap.h>
- #include <octave/gripes.h>
- #include <octave/lo-ieee.h>
--- 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 <config.h>
--#endif
--
- #include "base-lu.h"
- 
- template <class lu_type>
-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 <octave/config.h>
--#incude <octave/ov.h>
--])
--AC_CHECK_HEADER([octave/ls-oct-text.h], [PKG_CPPFLAGS="$PKG_CPPFLAGS -DHAVE_OCTAVE_TEXT_H=1"], [], 
--[
--#include <octave/config.h>
--#include <octave/ov.h>
--])
--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 <octave/pr-output.h>
- 
- 
--#if defined(HAVE_OCTAVE_TEXT_H)
--#  include <octave/ls-oct-text.h>
--#else
--#  include <octave/ls-oct-ascii.h>
--#endif
-+#include <octave/ls-oct-text.h>
- 
- #include "galois.h"
- #include "ov-galois.h"
--- 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
-+<http://www.gnu.org/licenses/>.
-+
-+*/
-+
-+#include "base-lu.h"
-+
-+template <class lu_type>
-+base_lu<lu_type>::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 <class lu_type>
-+bool
-+base_lu <lu_type> :: packed (void) const
-+{
-+  return l_fact.dims () == dim_vector ();
-+}
-+
-+template <class lu_type>
-+void
-+base_lu <lu_type> :: unpack (void)
-+{
-+  if (packed ())
-+    {
-+      l_fact = L ();
-+      a_fact = U (); // FIXME: sub-optimal
-+      ipvt = getp ();
-+    }
-+}
-+
-+template <class lu_type>
-+lu_type
-+base_lu <lu_type> :: 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 <class lu_type>
-+lu_type
-+base_lu <lu_type> :: 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 <class lu_type>
-+lu_type
-+base_lu <lu_type> :: Y (void) const
-+{
-+  if (! packed ())
-+    (*current_liboctave_error_handler)
-+      ("lu: Y () not implemented for unpacked form");
-+  return a_fact;
-+}
-+
-+template <class lu_type>
-+Array<octave_idx_type>
-+base_lu <lu_type> :: getp (void) const
-+{
-+  if (packed ())
-+    {
-+      octave_idx_type a_nr = a_fact.rows ();
-+
-+      Array<octave_idx_type> 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 <class lu_type>
-+PermMatrix
-+base_lu <lu_type> :: P (void) const
-+{
-+  return PermMatrix (getp (), false);
-+}
-+
-+template <class lu_type>
-+ColumnVector
-+base_lu <lu_type> :: P_vec (void) const
-+{
-+  octave_idx_type a_nr = a_fact.rows ();
-+
-+  ColumnVector p (a_nr);
-+
-+  Array<octave_idx_type> pvt = getp ();
-+
-+  for (octave_idx_type i = 0; i < a_nr; i++)
-+    p.xelem (i) = static_cast<double> (pvt.xelem (i) + 1);
-+
-+  return p;
-+}
-+
-+template <class lu_type>
-+bool
-+base_lu<lu_type>::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<int>& x, const int& n)
- {
- 
--  int x_len = x.length ();
-+  int x_len = x.numel ();
-   Array<int> si (dim_vector (n, 1), 0);
-   Array<int> 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
-+// <http://www.gnu.org/licenses/>.
-+//
-+// 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 <string>
-+
-+#include <octave/oct.h>
-+
-+// A simplified version of the filter function for specific lengths of a and b
-+// in the Galois field GF(2)
-+Array<int>
-+filter_gf2 (const Array<int>& b, const Array<int>& a,
-+            const Array<int>& x, const int& n)
-+{
-+
-+  int x_len = x.length ();
-+  Array<int> si (dim_vector (n, 1), 0);
-+  Array<int> 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<int>& a, const int& n, const int& m)
-+{
-+  Array<int> y (dim_vector (n+1, 1), 0);
-+  Array<int> x (dim_vector (n-m+2, 1), 0);
-+  y(0) = 1;
-+  y(n) = 1;
-+  x(0) = 1;
-+
-+  Array<int> b = filter_gf2 (y, a, x, n);
-+  b.resize (dim_vector (n+1, 1), 0);
-+  Array<int> p (dim_vector (m+1, 1), 0);
-+  p(0) = 1;
-+  Array<int> 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<int> 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<<i) ? 1 : 0);
-+    }
-+  else
-+    {
-+      Matrix tmp = args(1).matrix_value ();
-+      if ((tmp.rows () != 1) && (tmp.columns () != 1))
-+        {
-+          error ("cyclgen: generator polynomial must be a vector");
-+          return retval;
-+        }
-+
-+      if (tmp.rows () == 1)
-+        {
-+          mm = tmp.columns ();
-+          for (int j = 0; j < mm; j++) {
-+            if (tmp(0, j) == 1) {
-+              p |= ((unsigned long long)1 << j);
-+              pp(j) = 1;
-+            }
-+            else if (tmp(0, j) != 0) {
-+              error ("cyclgen: illegal generator polynomial");
-+              return retval;
-+            }
-+          }
-+        }
-+      else
-+        {
-+          mm = tmp.rows ();
-+          for (int i = 0; i < mm; i++)
-+            {
-+              if (tmp(i, 0) == 1)
-+                {
-+                  p |= ((unsigned long long)1 << i);
-+                  pp(i) = 1;
-+                }
-+              else if (tmp(i, 0) != 0)
-+                {
-+                  error ("cyclgen: illegal generator polynomial");
-+                  return retval;
-+                }
-+            }
-+        }
-+      mm = mm - 1;
-+    }
-+  k = n - mm;
-+
-+  if (nargin > 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<<mm))
-+        mask ^= p;
-+    }
-+
-+  Matrix parity (mm, n, 0);
-+  for (int i = 0; i < n; i++)
-+    for (int j = 0; j < mm; j++)
-+      if (alpha_to[i] & ((unsigned long long)1<<j))
-+        parity(j, i) = 1;
-+
-+  free (alpha_to);
-+  retval(0) = octave_value (parity);
-+
-+  if (nargout > 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<int>& x, const int& n)
- {
- 
--  int x_len = x.length ();
-+  int x_len = x.numel ();
-   Array<int> si (dim_vector (n, 1), 0);
-   Array<int> y (dim_vector (x_len, 1), 0);
- 
-@@ -217,8 +217,8 @@
-       for (unsigned long long i = (1UL<<m)+1; i < (1UL<<(1+m)); i+=2)
-         if (do_is_cyclic_polynomial (i, n, m))
-           {
--            cyclic_polys.resize (cyclic_polys.length ()+1);
--            cyclic_polys(cyclic_polys.length ()-1) = (double)i;
-+            cyclic_polys.resize (cyclic_polys.numel ()+1);
-+            cyclic_polys(cyclic_polys.numel ()-1) = (double)i;
-           }
-       break;
-     case CYCLIC_POLY_L:
-@@ -233,8 +233,8 @@
-             {
-               if (do_is_cyclic_polynomial (i, n, m))
-                 {
--                  cyclic_polys.resize (cyclic_polys.length ()+1);
--                  cyclic_polys(cyclic_polys.length ()-1) = (double)i;
-+                  cyclic_polys.resize (cyclic_polys.numel ()+1);
-+                  cyclic_polys(cyclic_polys.numel ()-1) = (double)i;
-                 }
-             }
-         }
-@@ -244,7 +244,7 @@
-       break;
-     }
- 
--  if (cyclic_polys.length () == 0)
-+  if (cyclic_polys.numel () == 0)
-     {
-       octave_stdout <<
-         "cyclpoly: no generator polynomial statifies constraints" << std::endl;
-@@ -254,8 +254,8 @@
-     {
-       if (polyrep)
-         {
--          Matrix polys (cyclic_polys.length (), m+1, 0);
--          for (int i = 0 ; i < cyclic_polys.length (); i++)
-+          Matrix polys (cyclic_polys.numel (), m+1, 0);
-+          for (int i = 0 ; i < cyclic_polys.numel (); i++)
-             for (int j = 0; j < m+1; j++)
-               if ((unsigned long long)cyclic_polys(i) & (1<<j))
-                 polys(i, j) = 1;
-diff -uNr a/src/cyclpoly.cc~ b/src/cyclpoly.cc~
---- a/src/cyclpoly.cc~	1969-12-31 19:00:00.000000000 -0500
-+++ b/src/cyclpoly.cc~	2018-04-09 13:52:05.442966638 -0400
-@@ -0,0 +1,282 @@
-+//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
-+// <http://www.gnu.org/licenses/>.
-+//
-+// 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 <iostream>
-+#include <string>
-+
-+#include <octave/oct.h>
-+
-+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<int>
-+filter_gf2 (const Array<int>& b, const Array<int>& a,
-+            const Array<int>& x, const int& n)
-+{
-+
-+  int x_len = x.numel ();
-+  Array<int> si (dim_vector (n, 1), 0);
-+  Array<int> 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<int> a (dim_vector (n+1, 1), 0);
-+  Array<int> y (dim_vector (n+1, 1), 0);
-+  Array<int> 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<int> b = filter_gf2 (y, a, x, n);
-+  b.resize(dim_vector (n+1, 1), 0);
-+  Array<int> p (dim_vector (m+1, 1), 0);
-+  p(0) = 1;
-+  Array<int> 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<<m)+1; i < (1UL<<(1+m)); i+=2)
-+        if (do_is_cyclic_polynomial (i, n, m))
-+          {
-+            cyclic_polys(0) = (double)i;
-+            break;
-+          }
-+      break;
-+    case CYCLIC_POLY_MAX:
-+      cyclic_polys.resize (1);
-+      for (unsigned long long i = (1UL<<(m+1))-1; i > (1UL<<m); i-=2)
-+        if (do_is_cyclic_polynomial (i, n, m))
-+          {
-+            cyclic_polys(0) = (double)i;
-+            break;
-+          }
-+      break;
-+    case CYCLIC_POLY_ALL:
-+      for (unsigned long long i = (1UL<<m)+1; i < (1UL<<(1+m)); i+=2)
-+        if (do_is_cyclic_polynomial (i, n, m))
-+          {
-+            cyclic_polys.resize (cyclic_polys.length ()+1);
-+            cyclic_polys(cyclic_polys.length ()-1) = (double)i;
-+          }
-+      break;
-+    case CYCLIC_POLY_L:
-+      for (unsigned long long i = ((unsigned long long)1<<m)+1;
-+           i < ((unsigned long long)1<<(1+m)); i+=2)
-+        {
-+          int li = 0;
-+          for (int j=0; j < m+1; j++)
-+            if (i & ((unsigned long long)1 << j))
-+              li++;
-+          if (li == l)
-+            {
-+              if (do_is_cyclic_polynomial (i, n, m))
-+                {
-+                  cyclic_polys.resize (cyclic_polys.length ()+1);
-+                  cyclic_polys(cyclic_polys.length ()-1) = (double)i;
-+                }
-+            }
-+        }
-+      break;
-+    default:
-+      error ("cyclpoly: impossible");
-+      break;
-+    }
-+
-+  if (cyclic_polys.length () == 0)
-+    {
-+      octave_stdout <<
-+        "cyclpoly: no generator polynomial statifies constraints" << std::endl;
-+      retval = octave_value (Matrix (0, 0));
-+    }
-+  else
-+    {
-+      if (polyrep)
-+        {
-+          Matrix polys (cyclic_polys.length (), m+1, 0);
-+          for (int i = 0 ; i < cyclic_polys.length (); i++)
-+            for (int j = 0; j < m+1; j++)
-+              if ((unsigned long long)cyclic_polys(i) & (1<<j))
-+                polys(i, j) = 1;
-+          retval = octave_value (polys);
-+        }
-+      else
-+        retval = octave_value (cyclic_polys);
-+    }
-+
-+  return retval;
-+}
-+
-+/*
-+%% Test input validation
-+%!error cyclpoly ()
-+%!error cyclpoly (1)
-+%!error cyclpoly (1, 2, 3, 4, 5)
-+*/
-+
-+/*
-+;;; Local Variables: ***
-+;;; mode: C++ ***
-+;;; End: ***
-+*/
-diff -uNr a/src/galois.cc b/src/galois.cc
---- a/src/galois.cc	2018-04-09 13:25:42.880981256 -0400
-+++ b/src/galois.cc	2018-04-09 13:53:27.547125644 -0400
-@@ -19,7 +19,7 @@
- // Initiative (www.opensource.org)
- 
- #include <octave/error.h>
--#include <octave/gripes.h>
-+#include <octave/errwarn.h>
- #include <octave/mx-op-defs.h>
- 
- #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<int> 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
-+// <http://www.gnu.org/licenses/>.
-+//
-+// 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 <octave/error.h>
-+#include <octave/errwarn.h>
-+#include <octave/mx-op-defs.h>
-+
-+#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<int>& a, const int& _m,
-+                const int& _primpoly) : MArray<int> (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<int>& a, const int& _m,
-+                const int& _primpoly) : MArray<int> (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<int> (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<int> (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<int> (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<int> (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<int>::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<int>::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<int>::index(i, j, resize_ok, rfv), m (), primpoly ());
-+
-+  return retval;
-+}
-+
-+galois
-+galois::concat (const galois& rb, const Array<int>& ra_idx)
-+{
-+  if (rb.numel () > 0)
-+    insert (rb, ra_idx(0), ra_idx(1));
-+  return *this;
-+}
-+
-+galois
-+galois::concat (const Matrix& rb, const Array<int>& ra_idx)
-+{
-+  if (numel () == 1)
-+    return *this;
-+
-+  galois tmp (0, 0, 0, m (), primpoly ());
-+  int _n = (1<<m ()) - 1;
-+  int r = rb.rows ();
-+  int c = rb.columns ();
-+  tmp.resize (dim_vector (r, c));
-+
-+  // Check the validity of the data in the matrix
-+  for (int i = 0; i < r; i++)
-+    {
-+      for (int j = 0; j < c; j++)
-+        {
-+          if ((rb(i, j) < 0) || (rb(i, j) > _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<int>& ra_idx)
-+{
-+  galois retval (0, 0, 0, rb.m (), rb.primpoly ());
-+  int _n = (1<<rb.m ()) - 1;
-+  int r = ra.rows ();
-+  int c = ra.columns ();
-+  retval.resize (dim_vector (r, c));
-+  if (ra.numel () < 1)
-+    return retval;
-+
-+  // FIXME:
-+  // Check the validity of the data in the matrix. This is problematic
-+  // as "ra" is not initialized on the initial resize and so contains
-+  // random data that will be replaced. Humm, disable for now
-+  for (int i = 0; i < r; i++)
-+    {
-+      for (int j = 0; j < c; j++)
-+        {
-+#if 0
-+          if ((ra(i, j) < 0) || (ra(i, j) > _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<int>::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<bool, int> (*this, dim, mx_inline_all);
-+}
-+
-+boolMatrix
-+galois::any (int dim) const
-+{
-+  return do_mx_red_op<bool, int> (*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 <galois>;
-+
-+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<int> 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<int> 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 <octave/defun-dld.h>
--#include <octave/gripes.h>
-+#include <octave/errwarn.h>
-+#include <octave/interpreter.h>
- #include <octave/oct-locbuf.h>
- #include <octave/ov.h>
- #include <octave/utils.h>
-@@ -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 <owner@eccpage.com>
-+// Copyright (C) 2002 Phil Karn <karn@ka9q.net>
-+// 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
-+// <http://www.gnu.org/licenses/>.
-+//
-+// 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 <octave/defun-dld.h>
-+#include <octave/errwarn.h>
-+#include <octave/interpreter.h>
-+#include <octave/oct-locbuf.h>
-+#include <octave/ov.h>
-+#include <octave/utils.h>
-+#include <octave/variables.h>
-+
-+#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<<m))
-+    m++;
-+  int nn = (1<<m) - 1;
-+
-+  if (msg.cols () != k)
-+    {
-+      error ("rsenc: message contains incorrect number of symbols");
-+      return retval;
-+    }
-+
-+  if (msg.m () != m)
-+    {
-+      error ("rsenc: message in incorrect galois field for codeword length");
-+      return retval;
-+    }
-+
-+  if ((n < 3) || (n < k) || (m > __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++)
-+            msg(l, i) = 0;
-+        }
-+      for (int l = 0; l < nsym; l++)
-+        {
-+          galois par (nroots, 1, 0, m, primpoly);
-+          for (int i = n; 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; i<nroots; i++)
-+          if(s[i] == 0)
-+            s[i] = data(drow, j);
-+          else
-+            s[i] = data(drow, j) ^ data.alpha_to (modn (data.index_of (s[i]) +
-+                                                        (fcr+i)*prim, m, n));
-+    }
-+  else
-+    {
-+      for (i = 0; i<nroots; i++)
-+        s[i] = data(drow, n-1);
-+
-+      for (j = n-1; j>0; 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(&reg[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<<m))
-+    m++;
-+  int nn = (1<<m) - 1;
-+
-+  if (code.cols () != n)
-+    {
-+      error ("rsdec: coded message contains incorrect number of symbols");
-+      return retval;
-+    }
-+
-+  if (code.m () != m)
-+    {
-+      error ("rsdec: coded message in incorrect galois field for "
-+             "codeword length");
-+      return retval;
-+    }
-+
-+  if ((n < 3) || (n < k) || (m > __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<nroots; i++)
-+                {
-+                  if (val == 0)
-+                    val = genpoly(i+1, 0);
-+                  else
-+                    val = genpoly(i+1, 0) ^ genpoly.alpha_to (modn (indx +
-+                                                                    genpoly.index_of (val),
-+                                                                    m, nn));
-+                }
-+              if (val == 0)
-+                {
-+                  roots[count] = j;
-+                  count++;
-+                  if (count == nroots)
-+                    break;
-+                }
-+            }
-+
-+          if (count != nroots)
-+            {
-+              error ("rsdec: generator polynomial can not have repeated roots");
-+              return retval;
-+            }
-+
-+          // Logarithm of roots wrt primitive element
-+          for (int i = 0; i < count; i++)
-+            roots[i] = genpoly.index_of (roots[i]);
-+
-+          // Find a corresponding fcr and prim that coincide with the roots.
-+          // FIXME: This is a naive algorithm and should be improved !!!
-+          bool found = true;
-+          for (fcr = 1; fcr < n+1; fcr++)
-+            {
-+              for (prim = 1; prim < n+1; prim++)
-+                {
-+                  found = true;
-+                  for (int i = 0; i<nroots; i++)
-+                    {
-+                      int tmp = modn ((fcr + i)*prim, m, n);
-+                      for (int j = 0; j<count; j++)
-+                        {
-+                          if (tmp == roots[j])
-+                            {
-+                              tmp = -1;
-+                              break;
-+                            }
-+                        }
-+                      if (tmp != -1)
-+                        {
-+                          found = false;
-+                          break;
-+                        }
-+                    }
-+                  if (found)
-+                    break;
-+                }
-+              if (found)
-+                break;
-+            }
-+        }
-+    }
-+  else
-+    {
-+      fcr = 1;
-+      prim = 1;
-+    }
-+
-+  /* Find prim-th root of 1, used in decoding */
-+  for (iprim = 1; (iprim % prim) != 0; iprim += n)
-+    ;
-+  iprim = iprim / prim;
-+
-+  galois msg (nsym, k, 0, m, primpoly);
-+  ColumnVector nerr (nsym, 0);
-+
-+  if (nn != n)
-+    {
-+      code.resize (dim_vector (nsym, nn), 0);
-+      if (parity_at_end)
-+        for (int l = 0; l < nsym; l++)
-+          for (int i=n; 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<<m))
-+    m++;
-+
-+  int n = (1<<m) - 1;
-+
-+  if (msg.cols () != k)
-+    {
-+      error ("bchenco: message contains incorrect number of symbols");
-+      return retval;
-+    }
-+
-+  if ((n < 3) || (nn < k) || (m > __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<int> 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<n; i++)
-+            if ((found(i) == 0) && (c.index_of (i+1) < idx))
-+              idx = c.index_of (i+1);
-+
-+          c.resize (dim_vector (nc+1, m));
-+          cs.resize (dim_vector (nc+1, 1));
-+          c(nc, 0) = idx;
-+          found(c.alpha_to (idx)-1) = 1;
-+          cs(nc) = 1;
-+          int r = idx;
-+          while ((r = modn (r<<1, m, n)) > 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<nn-k; i++)
-+            msg(l, i) = 0;
-+        }
-+
-+      for (int l = 0; l < nsym; l++)
-+        {
-+          for (int i = k-1; 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<<m))
-+    m++;
-+
-+  int n = (1<<m) - 1;
-+
-+  if ((n < 3) || (n < k) || (m > __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<<i);
-+            }
-+        }
-+    }
-+
-+  // Create a variable in the require Galois Field to have access to the
-+  // lookup tables alpha_to and index_of.
-+  galois tables (1, 1, 0, m, prim);
-+  ColumnVector nerr (nsym, 0);
-+
-+  for (int lsym = 0; lsym < nsym; lsym++)
-+    {
-+      /* first form the syndromes */
-+      Array<int> 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<int> 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<int> 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<char> 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
-+// <http://www.gnu.org/licenses/>.
-+//
-+// 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 <octave/oct.h>
-+
-+static int
-+get_weight (const Array<char>& 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<char> 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<char> 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 <iostream>
- 
- #include <octave/byte-swap.h>
--#include <octave/gripes.h>
-+#include <octave/errwarn.h>
- #include <octave/lo-ieee.h>
- #include <octave/oct-locbuf.h>
- #include <octave/oct-obj.h>
-@@ -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
-+// <http://www.gnu.org/licenses/>.
-+//
-+// 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 <octave/ov.h>
-+#include <octave/ov-typeinfo.h>
-+
-+#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<octave_value_list>& idx);
-+
-+  octave_value_list subsref (const std::string& type,
-+                             const std::list<octave_value_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: ***
-+*/
--- 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 <octave/mx-base.h>
-+#include <octave/ov-typeinfo.h>
- 
- #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 <octave/oct-obj.h>
-+#include <octave/ovl.h>
- #include <octave/ops.h>
- 
- #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 <octave/errwarn.h>
- #include <octave/lo-ieee.h>
- #include <octave/oct-locbuf.h>
--#include <octave/oct-obj.h>
- #include <octave/ov.h>
-+#include <octave/ovl.h>
- #include <octave/pr-output.h>
- 
- 
-@@ -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<<m)+1; i < (1<<(1+m)); i+=2)
-         if (do_isprimitive (i, m))
-           {
--            primpolys.resize (primpolys.length ()+1);
--            primpolys(primpolys.length ()-1) = (double)i;
-+            primpolys.resize (primpolys.numel ()+1);
-+            primpolys(primpolys.numel ()-1) = (double)i;
-           }
-       break;
-     case PRIMPOLY_K:
-@@ -231,8 +231,8 @@
-             {
-               if (do_isprimitive (i, m))
-                 {
--                  primpolys.resize (primpolys.length ()+1);
--                  primpolys(primpolys.length ()-1) = (double)i;
-+                  primpolys.resize (primpolys.numel ()+1);
-+                  primpolys(primpolys.numel ()-1) = (double)i;
-                 }
-             }
-         }
-@@ -244,12 +244,12 @@
- 
-   if (disp)
-     {
--      if (primpolys.length () == 0)
-+      if (primpolys.numel () == 0)
-         warning ("primpoly: No primitive polynomial satisfies the given constraints");
-       else
-         {
-           octave_stdout << std::endl << "Primitive polynomial(s) =" << std::endl << std::endl;
--          for (int i = 0; i < primpolys.length (); i++)
-+          for (int i = 0; i < primpolys.numel (); i++)
-             {
-               bool first = true;
-               for (int j = m; j >= 0; j--)
-diff -r old src/ov-galois.cc
---- a/src/ov-galois.cc	2019-06-06 17:27:40.000000000 -0400
-+++ b/src/ov-galois.cc	2019-06-07 07:18:23.978697671 -0400
-@@ -21,6 +21,7 @@
- #include <iostream>
- 
- #include <octave/byte-swap.h>
-+#include <octave/error.h>
- #include <octave/errwarn.h>
- #include <octave/lo-ieee.h>
- #include <octave/oct-locbuf.h>
--- a/src/of-communications.mk	Wed Jan 01 18:52:37 2020 -0500
+++ b/src/of-communications.mk	Mon Jan 06 06:19:36 2020 -0500
@@ -3,8 +3,8 @@
 
 PKG             := of-communications
 $(PKG)_IGNORE   :=
-$(PKG)_VERSION  := 1.2.1
-$(PKG)_CHECKSUM := bf70d8c315c2239e168c02522482c81e4b912968
+$(PKG)_VERSION  := 1.2.2
+$(PKG)_CHECKSUM := 90ebf5cd84ba8df1f8c14241598d0baedc39f371
 $(PKG)_REMOTE_SUBDIR := 
 $(PKG)_SUBDIR   := communications-$($(PKG)_VERSION)
 $(PKG)_FILE     := communications-$($(PKG)_VERSION).tar.gz