Mercurial > hg > octave-max

changeset 237:5a9e23307fb0

Find changesets by keywords (author, files, the commit message), revision number or hash, or revset expression.
[project @ 1993-11-30 20:23:04 by jwe] Initial revision
author jwe
date Tue, 30 Nov 1993 20:23:04 +0000
parents 93bbd907de34
children 780cbbc57b7c
files liboctave/Array.cc liboctave/MArray.cc liboctave/MArray.h liboctave/mx-kludge.cc liboctave/mx-kludge.h
diffstat 5 files changed, 2336 insertions(+), 0 deletions(-) [+]
line wrap: on
line diff
new file mode 100644
--- /dev/null
+++ b/liboctave/Array.cc
@@ -0,0 +1,938 @@
+// Template array classes                              -*- C++ -*-
+/*
+
+Copyright (C) 1993 John W. Eaton
+
+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 2, 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, write to the Free
+Software Foundation, 675 Mass Ave, Cambridge, MA 02139, USA.
+
+*/
+
+#ifdef HAVE_CONFIG_H
+#include "config.h"
+#endif
+
+#include <assert.h>
+
+#if defined (__GNUG__) && defined (USE_EXTERNAL_TEMPLATES)
+#pragma implementation
+#endif
+
+#include "Array.h"
+#include "lo-error.h"
+
+/*
+ * The real representation of all arrays.
+ */
+
+template <class T>
+ArrayRep<T>::ArrayRep (T *d, int l)
+{
+  data = d;
+  len = l;
+}
+
+template <class T>
+ArrayRep<T>::ArrayRep (void)
+{
+  len = 0;
+  data = (T *) 0;
+}
+
+template <class T>
+ArrayRep<T>::ArrayRep (int n)
+{
+  len = n;
+  data = new T [len];
+}
+
+template <class T>
+ArrayRep<T>::ArrayRep (const ArrayRep<T>& a)
+{
+  len = a.len;
+  count = a.count;
+  data = new T [len];
+  for (int i = 0; i < len; i++)
+    data[i] = a.data[i];
+}
+
+template <class T>
+ArrayRep<T>::~ArrayRep (void)
+{
+  delete [] data;
+  data = (T *) 0;
+}
+
+template <class T>
+int
+ArrayRep<T>::length (void) const
+{
+  return len;
+}
+
+template <class T>
+T&
+ArrayRep<T>::elem (int n)
+{
+  return data[n];
+}
+
+template <class T>
+T
+ArrayRep<T>::elem (int n) const
+{
+  return data[n];
+}
+
+/*
+ * One dimensional array class.  Handles the reference counting for
+ * all the derived classes.
+ */
+
+template <class T>
+Array<T>::Array (T *d, int l)
+{
+  rep = new ArrayRep<T> (d, l);
+  rep->count = 1;
+}
+
+template <class T>
+Array<T>::Array (void)
+{
+  rep = new ArrayRep<T>;
+  rep->count = 1;
+}
+
+template <class T>
+Array<T>::Array (int n)
+{
+  rep = new ArrayRep<T> (n);
+  rep->count = 1;
+}
+
+template <class T>
+Array<T>::Array (int n, const T& val)
+{
+  rep = new ArrayRep<T> (n);
+  rep->count = 1;
+  for (int i = 0; i < n; i++)
+    rep->data[i] = val;
+}
+
+template <class T>
+Array<T>::Array (const Array<T>& a)
+{
+  rep = a.rep;
+  rep->count++;
+}
+
+template <class T>
+Array<T>::~Array (void)
+{
+  if (--rep->count <= 0)
+    delete rep;
+}
+
+template <class T>
+Array<T>&
+Array<T>::operator = (const Array<T>& a)
+{
+  if (--rep->count <= 0)
+    delete rep;
+
+  rep = a.rep;
+  rep->count++;
+  return *this;
+}
+
+template <class T>
+int
+Array<T>::capacity (void) const
+{
+  return rep->length ();
+}
+
+template <class T>
+int
+Array<T>::length (void) const
+{
+  return rep->length ();
+}
+
+template <class T>
+T&
+Array<T>::elem (int n)
+{
+  if (rep->count > 1)
+    {
+      --rep->count;
+      rep = new ArrayRep<T> (*rep);
+      rep->count = 1;
+    }
+  return rep->elem (n);
+}
+
+template <class T>
+T&
+Array<T>::checkelem (int n)
+{
+  if (n < 0 || n >= rep->length ())
+    {
+      (*current_liboctave_error_handler) ("range error");
+      static T foo (0);
+      return foo;
+    }
+  return elem (n);
+}
+
+template <class T>
+T&
+Array<T>::operator () (int n)
+{
+  return checkelem (n);
+}
+
+template <class T>
+T&
+Array<T>::xelem (int n)
+{
+  return rep->elem (n);
+}
+
+template <class T>
+T
+Array<T>::elem (int n) const
+{
+  return rep->elem (n);
+}
+
+template <class T>
+T
+Array<T>::checkelem (int n) const
+{
+  if (n < 0 || n >= rep->length ())
+    {
+      (*current_liboctave_error_handler) ("range error");
+      return T (0);
+    }
+  return elem (n);
+}
+
+template <class T>
+T
+Array<T>::operator () (int n) const
+{
+  return checkelem (n);
+}
+
+template <class T>
+void
+Array<T>::resize (int n)
+{
+  if (n < 0)
+    {
+      (*current_liboctave_error_handler)
+	("can't resize to negative dimension");
+      return;
+    }
+
+  if (n == length ())
+    return;
+
+  ArrayRep<T> *old_rep = rep;
+  const T *old_data = data ();
+  int old_len = length ();
+
+  rep = new ArrayRep<T> (n);
+  rep->count = 1;
+
+  if (old_data && old_len > 0)
+    {
+      int min_len = old_len < n ? old_len : n;
+
+      for (int i = 0; i < min_len; i++)
+	xelem (i) = old_data[i];
+    }
+
+  if (--old_rep->count <= 0)
+    delete old_rep;
+}
+
+template <class T>
+void
+Array<T>::resize (int n, const T& val)
+{
+  if (n < 0)
+    {
+      (*current_liboctave_error_handler)
+	("can't resize to negative dimension");
+      return;
+    }
+
+  if (n == length ())
+    return;
+
+  ArrayRep<T> *old_rep = rep;
+  const T *old_data = data ();
+  int old_len = length ();
+
+  rep = new ArrayRep<T> (n);
+  rep->count = 1;
+
+  int min_len = old_len < n ? old_len : n;
+
+  if (old_data && old_len > 0)
+    {
+      for (int i = 0; i < min_len; i++)
+	xelem (i) = old_data[i];
+    }
+
+  for (int i = old_len; i < n; i++)
+    xelem (i) = val;
+
+  if (--old_rep->count <= 0)
+    delete old_rep;
+}
+
+template <class T>
+const T *
+Array<T>::data (void) const
+{
+  return rep->data;
+}
+
+template <class T>
+T *
+Array<T>::fortran_vec (void)
+{
+  if (rep->count > 1)
+    {
+      --rep->count;
+      rep = new ArrayRep<T> (*rep);
+      rep->count = 1;
+    }
+  return rep->data;
+}
+
+/*
+ * Two dimensional array class.
+ */
+
+template <class T>
+Array2<T>::Array2 (T *d, int n, int m) : Array<T> (d, n*m)
+{
+  d1 = n;
+  d2 = m;
+}
+
+template <class T>
+Array2<T>::Array2 (void) : Array<T> ()
+{
+  d1 = 0;
+  d2 = 0;
+}
+
+template <class T>
+Array2<T>::Array2 (int n, int m) : Array<T> (n*m)
+{
+  d1 = n;
+  d2 = m;
+}
+
+template <class T>
+Array2<T>::Array2 (int n, int m, const T& val) : Array<T> (n*m, val)
+{
+  d1 = n;
+  d2 = m;
+}
+
+template <class T>
+Array2<T>::Array2 (const Array2<T>& a) : Array<T> (a)
+{
+  d1 = a.d1;
+  d2 = a.d2;
+}
+
+template <class T>
+Array2<T>::Array2 (const DiagArray<T>& a)
+  : Array<T> (a.rows () * a.cols (), T (0))
+{
+  for (int i = 0; i < a.length (); i++)
+    elem (i, i) = a.elem (i, i);
+}
+
+template <class T>
+Array2<T>&
+Array2<T>::operator = (const Array2<T>& a)
+{
+  Array<T>::operator = (a);
+  d1 = a.d1;
+  d2 = a.d2;
+  return *this;
+}
+
+template <class T>
+int
+Array2<T>::dim1 (void) const
+{
+  return d1;
+}
+
+template <class T>
+int
+Array2<T>::dim2 (void) const
+{
+  return d2;
+}
+
+template <class T>
+int
+Array2<T>::rows (void) const
+{
+  return d1;
+}
+
+template <class T>
+int
+Array2<T>::cols (void) const
+{
+  return d2;
+}
+
+template <class T>
+int
+Array2<T>::columns (void) const
+{
+  return d2;
+}
+
+template <class T>
+T&
+Array2<T>::elem (int i, int j)
+{
+  return Array<T>::elem (d1*j+i);
+}
+
+template <class T>
+T&
+Array2<T>::checkelem (int i, int j)
+{
+  return Array<T>::checkelem (d1*j+i);
+}
+
+template <class T>
+T&
+Array2<T>::operator () (int i, int j)
+{
+  return Array<T>::checkelem (d1*j+i);
+}
+
+template <class T>
+T&
+Array2<T>::xelem (int i, int j)
+{
+  return Array<T>::xelem (d1*j+i);
+}
+
+template <class T>
+T
+Array2<T>::elem (int i, int j) const
+{
+  return Array<T>::elem (d1*j+i);
+}
+
+template <class T>
+T
+Array2<T>::checkelem (int i, int j) const
+{
+  return Array<T>::checkelem (d1*j+i);
+}
+
+template <class T>
+T
+Array2<T>::operator () (int i, int j) const
+{
+  return Array<T>::checkelem (d1*j+i);
+}
+
+template <class T>
+void
+Array2<T>::resize (int r, int c)
+{
+  if (r < 0 || c < 0)
+    {
+      (*current_liboctave_error_handler)
+	("can't resize to negative dimension");
+      return;
+    }
+
+  if (r == dim1 () && c == dim2 ())
+    return;
+
+  ArrayRep<T> *old_rep = rep;
+  const T *old_data = data ();
+  int old_d1 = dim1 ();
+  int old_d2 = dim2 ();
+  int old_len = length ();
+
+  rep = new ArrayRep<T> (r*c);
+  rep->count = 1;
+
+  d1 = r;
+  d2 = c;
+
+  if (old_data && old_len > 0)
+    {
+      int min_r = old_d1 < r ? old_d1 : r;
+      int min_c = old_d2 < c ? old_d2 : c;
+
+      for (int j = 0; j < min_c; j++)
+	for (int i = 0; i < min_r; i++)
+	  xelem (i, j) = old_data[old_d1*j+i];
+    }
+
+  if (--old_rep->count <= 0)
+    delete old_rep;
+}
+
+template <class T>
+void
+Array2<T>::resize (int r, int c, const T& val)
+{
+  if (r < 0 || c < 0)
+    {
+      (*current_liboctave_error_handler)
+	("can't resize to negative dimension");
+      return;
+    }
+
+  if (r == dim1 () && c == dim2 ())
+    return;
+
+  ArrayRep<T> *old_rep = rep;
+  const T *old_data = data ();
+  int old_d1 = dim1 ();
+  int old_d2 = dim2 ();
+  int old_len = length ();
+
+  rep = new ArrayRep<T> (r*c);
+  rep->count = 1;
+
+  d1 = r;
+  d2 = c;
+
+  int min_r = old_d1 < r ? old_d1 : r;
+  int min_c = old_d2 < c ? old_d2 : c;
+
+  int i, j;
+
+  if (old_data && old_len > 0)
+    {
+      for (j = 0; j < min_c; j++)
+	for (i = 0; i < min_r; i++)
+	  xelem (i, j) = old_data[old_d1*j+i];
+    }
+
+  for (j = 0; j < min_c; j++)
+    for (i = min_r; i < r; i++)
+      xelem (i, j) = val;
+
+  for (j = min_c; j < c; j++)
+    for (i = 0; i < r; i++)
+      xelem (i, j) = val;
+
+  if (--old_rep->count <= 0)
+    delete old_rep;
+}
+
+/*
+ * Three dimensional array class.
+ */
+
+template <class T>
+Array3<T>::Array3 (T *d, int n, int m, int k) : Array2<T> (d, n, m*k)
+{
+  d2 = m;
+  d3 = k;
+}
+
+template <class T>
+Array3<T>::Array3 (void) : Array2<T> ()
+{
+  d2 = 0;
+  d3 = 0;
+}
+
+template <class T>
+Array3<T>::Array3 (int n, int m, int k) : Array2<T> (n, m*k)
+{
+  d2 = m;
+  d3 = k;
+}
+
+template <class T>
+Array3<T>::Array3 (int n, int m, int k, const T& val) : Array2<T> (n, m*k, val)
+{
+  d2 = m;
+  d3 = k;
+}
+
+template <class T>
+Array3<T>::Array3 (const Array3<T>& a) : Array2<T> (a)
+{
+  d2 = a.d2;
+  d3 = a.d3;
+}
+
+template <class T>
+Array3<T>&
+Array3<T>::operator = (const Array3<T>& a)
+{
+  Array<T>::operator = (a);
+  d1 = a.d1;
+  d2 = a.d2;
+  d3 = a.d3;
+  return *this;
+}
+
+template <class T>
+int
+Array3<T>::dim3 (void) const
+{
+  return d3;
+}
+
+template <class T>
+T&
+Array3<T>::elem (int i, int j, int k)
+{
+  return Array2<T>::elem (i, d2*k+j);
+}
+
+template <class T>
+T&
+Array3<T>::checkelem (int i, int j, int k)
+{
+  return Array2<T>::checkelem (i, d1*k+j);
+}
+
+template <class T>
+T&
+Array3<T>::operator () (int i, int j, int k)
+{
+  return Array2<T>::checkelem (i, d2*k+j);
+}
+
+template <class T>
+T&
+Array3<T>::xelem (int i, int j, int k)
+{
+  return Array2<T>::xelem (i, d2*k+j);
+}
+
+template <class T>
+T
+Array3<T>::elem (int i, int j, int k) const
+{
+  return Array2<T>::elem (i, d2*k+j);
+}
+
+template <class T>
+T
+Array3<T>::checkelem (int i, int j, int k) const
+{
+  return Array2<T>::checkelem (i, d1*k+j);
+}
+
+template <class T>
+T
+Array3<T>::operator () (int i, int j, int k) const
+{
+  return Array2<T>::checkelem (i, d2*k+j);
+}
+
+template <class T>
+void
+Array3<T>::resize (int n, int m, int k)
+{
+  assert (0); /* XXX FIXME XXX */
+}
+
+template <class T>
+void
+Array3<T>::resize (int n, int m, int k, const T& val)
+{
+  assert (0); /* XXX FIXME XXX */
+}
+
+/*
+ * A two-dimensional array with diagonal elements only.
+ */
+
+template <class T>
+DiagArray<T>::DiagArray (T *d, int r, int c) : Array<T> (d, r < c ? r : c)
+{
+  nr = r;
+  nc = c;
+}
+
+template <class T>
+DiagArray<T>::DiagArray (void) : Array<T> ()
+{
+  nr = 0;
+  nc = 0;
+}
+
+template <class T>
+DiagArray<T>::DiagArray (int n) : Array<T> (n)
+{
+  nr = n;
+  nc = n;
+}
+
+template <class T>
+DiagArray<T>::DiagArray (int n, const T& val) : Array<T> (n, val)
+{
+  nr = 0;
+  nc = 0;
+}
+
+template <class T>
+DiagArray<T>::DiagArray (int r, int c) : Array<T> (r < c ? r : c)
+{
+  nr = r;
+  nc = c;
+}
+
+template <class T>
+DiagArray<T>::DiagArray (int r, int c, const T& val)
+  : Array<T> (r < c ? r : c, val)
+{
+  nr = r;
+  nc = c;
+}
+
+template <class T>
+DiagArray<T>::DiagArray (const Array<T>& a) : Array<T> (a)
+{
+  nr = nc = a.length ();
+}
+
+template <class T>
+DiagArray<T>::DiagArray (const DiagArray<T>& a) : Array<T> (a)
+{
+  nr = a.nr;
+  nc = a.nc;
+}
+
+template <class T>
+DiagArray<T>&
+DiagArray<T>::operator = (const DiagArray<T>& a)
+{
+  Array<T>::operator = (a);
+  nr = a.nr;
+  nc = a.nc;
+  return *this;
+}
+
+template <class T>
+int
+DiagArray<T>::dim1 (void) const
+{
+  return nr;
+}
+
+template <class T>
+int
+DiagArray<T>::dim2 (void) const
+{
+  return nc;
+}
+
+template <class T>
+int
+DiagArray<T>::rows (void) const
+{
+  return nr;
+}
+
+template <class T>
+int
+DiagArray<T>::cols (void) const
+{
+  return nc;
+}
+
+template <class T>
+int
+DiagArray<T>::columns (void) const
+{
+  return nc;
+}
+
+template <class T>
+T&
+DiagArray<T>::elem (int r, int c) 
+{
+  static T foo (0);
+  return (r == c) ? Array<T>::elem (r) : foo;
+}
+
+template <class T>
+T&
+DiagArray<T>::checkelem (int r, int c) 
+{
+  static T foo (0);
+  return (r == c) ? Array<T>::checkelem (r) : foo;
+}
+
+template <class T>
+T&
+DiagArray<T>::operator () (int r, int c) 
+{
+  static T foo (0);
+  return (r == c) ? Array<T>::operator () (r) : foo;
+}
+
+template <class T>
+T&
+DiagArray<T>::xelem (int r, int c)
+{
+  static T foo (0);
+  return (r == c) ? Array<T>::xelem (r) : foo;
+}
+
+template <class T>
+T
+DiagArray<T>::elem (int r, int c) const
+{
+  return (r == c) ? Array<T>::elem (r) : T (0);
+}
+
+template <class T>
+T
+DiagArray<T>::checkelem (int r, int c) const
+{
+  return (r == c) ? Array<T>::checkelem (r) : T (0);
+}
+
+template <class T>
+T
+DiagArray<T>::operator () (int r, int c) const
+{
+  return (r == c) ? Array<T>::operator () (r) : T (0);
+}
+
+template <class T>
+void
+DiagArray<T>::resize (int r, int c)
+{
+  if (r < 0 || c < 0)
+    {
+      (*current_liboctave_error_handler)
+	("can't resize to negative dimensions");
+      return;
+    }
+
+  if (r == dim1 () && c == dim2 ())
+    return;
+
+  ArrayRep<T> *old_rep = rep;
+  const T *old_data = data ();
+  int old_len = length ();
+
+  int new_len = r < c ? r : c;
+
+  rep = new ArrayRep<T> (new_len);
+  rep->count = 1;
+
+  nr = r;
+  nc = c;
+
+  if (old_data && old_len > 0)
+    {
+      int min_len = old_len < new_len ? old_len : new_len;
+
+      for (int i = 0; i < min_len; i++)
+	xelem (i, i) = old_data[i];
+    }
+
+  if (--old_rep->count <= 0)
+    delete old_rep;
+}
+
+template <class T>
+void
+DiagArray<T>::resize (int r, int c, const T& val)
+{
+  if (r < 0 || c < 0)
+    {
+      (*current_liboctave_error_handler)
+	("can't resize to negative dimensions");
+      return;
+    }
+
+  if (r == dim1 () && c == dim2 ())
+    return;
+
+  ArrayRep<T> *old_rep = rep;
+  const T *old_data = data ();
+  int old_len = length ();
+
+  int new_len = r < c ? r : c;
+
+  rep = new ArrayRep<T> (new_len);
+  rep->count = 1;
+
+  nr = r;
+  nc = c;
+
+  int min_len = old_len < new_len ? old_len : new_len;
+
+  if (old_data && old_len > 0)
+    {
+      for (int i = 0; i < min_len; i++)
+	xelem (i, i) = old_data[i];
+    }
+
+  for (int i = min_len; i < new_len; i++)
+    xelem (i, i) = val;
+
+  if (--old_rep->count <= 0)
+    delete old_rep;
+}
+
+#ifdef __GNUG__
+#if defined (OCTAVE_SOURCE) && defined (USE_EXTERNAL_TEMPLATES)
+
+typedef Array<double> array_type_double;
+typedef Array2<double> array2_type_double;
+typedef DiagArray<double> diag_array_type_double;
+
+#include <Complex.h>
+typedef Array<Complex> array_type_complex;
+typedef Array2<Complex> array2_type_complex;
+typedef DiagArray<Complex> diag_array_type_complex;
+
+#endif
+#endif
+
+/*
+;;; Local Variables: ***
+;;; mode: C++ ***
+;;; page-delimiter: "^/\\*" ***
+;;; End: ***
+*/
new file mode 100644
--- /dev/null
+++ b/liboctave/MArray.cc
@@ -0,0 +1,545 @@
+// MArray.cc                                             -*- C++ -*-
+/*
+
+Copyright (C) 1992, 1993 John W. Eaton
+
+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 2, 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, write to the Free
+Software Foundation, 675 Mass Ave, Cambridge, MA 02139, USA.
+
+*/
+
+#ifdef HAVE_CONFIG_H
+#include "config.h"
+#endif
+
+#include "MArray.h"
+#include "lo-error.h"
+
+// Nothing like a little CPP abuse to brighten everyone's day.  Would
+// have been nice to do this with template functions but as of 2.5.x,
+// g++ seems to fail to resolve them properly.
+
+#define DO_VS_OP(OP) \
+  int l = a.length (); \
+  T *result = 0; \
+  if (l > 0) \
+    { \
+      result = new T [l]; \
+      const T *x = a.data (); \
+      for (int i = 0; i < l; i++) \
+	result[i] = x[i] OP s; \
+    }
+
+#define DO_SV_OP(OP) \
+  int l = a.length (); \
+  T *result = 0; \
+  if (l > 0) \
+    { \
+      result = new T [l]; \
+      const T *x = a.data (); \
+      for (int i = 0; i < l; i++) \
+	result[i] = s OP x[i]; \
+    }
+
+#define DO_VV_OP(OP) \
+  T *result = 0; \
+  if (l > 0) \
+    { \
+      result = new T [l]; \
+      const T *x = a.data (); \
+      const T *y = b.data (); \
+      for (int i = 0; i < l; i++) \
+	result[i] = x[i] OP y[i]; \
+    }
+
+#define NEG_V \
+  int l = a.length (); \
+  T *result = 0; \
+  if (l > 0) \
+    { \
+      result = new T [l]; \
+      const T *x = a.data (); \
+      for (int i = 0; i < l; i++) \
+	result[i] = -x[i]; \
+    }
+
+/*
+ * One dimensional array with math ops.
+ */
+
+// Element by element MArray by scalar ops.
+
+template <class T>
+MArray<T>
+operator + (const MArray<T>& a, const T& s)
+{
+  DO_VS_OP (+);
+  return MArray<T> (result, l);
+}
+
+template <class T>
+MArray<T>
+operator - (const MArray<T>& a, const T& s)
+{
+  DO_VS_OP (-);
+  return MArray<T> (result, l);
+}
+
+template <class T>
+MArray<T>
+operator * (const MArray<T>& a, const T& s)
+{
+  DO_VS_OP (*);
+  return MArray<T> (result, l);
+}
+
+template <class T>
+MArray<T>
+operator / (const MArray<T>& a, const T& s)
+{
+  DO_VS_OP (/);
+  return MArray<T> (result, l);
+}
+
+// Element by element scalar by MArray ops.
+
+template <class T>
+MArray<T>
+operator + (const T& s, const MArray<T>& a)
+{
+  DO_SV_OP (+);
+  return MArray<T> (result, l);
+}
+
+template <class T>
+MArray<T>
+operator - (const T& s, const MArray<T>& a)
+{
+  DO_SV_OP (-);
+  return MArray<T> (result, l);
+}
+
+template <class T>
+MArray<T>
+operator * (const T& s, const MArray<T>& a)
+{
+  DO_SV_OP (*);
+  return MArray<T> (result, l);
+}
+
+template <class T>
+MArray<T>
+operator / (const T& s, const MArray<T>& a)
+{
+  DO_SV_OP (/);
+  return MArray<T> (result, l);
+}
+
+// Element by element MArray by MArray ops.
+
+template <class T>
+MArray<T>
+operator + (const MArray<T>& a, const MArray<T>& b)
+{
+  int l = a.length ();
+  if (l != b.length ())
+    {
+      (*current_liboctave_error_handler)
+	("nonconformant array addition attempted");
+      return MArray<T> ();
+    }
+
+  if (l == 0)
+    return MArray<T> ();
+
+  DO_VV_OP (+);
+  return MArray<T> (result, l);
+}
+
+template <class T>
+MArray<T>
+operator - (const MArray<T>& a, const MArray<T>& b)
+{
+  int l = a.length ();
+  if (l != b.length ())
+    {
+      (*current_liboctave_error_handler)
+	("nonconformant array subtraction attempted");
+      return MArray<T> ();
+    }
+
+  if (l == 0)
+    return MArray<T> ();
+
+  DO_VV_OP (-);
+  return MArray<T> (result, l);
+}
+
+template <class T>
+MArray<T>
+product (const MArray<T>& a, const MArray<T>& b)
+{
+  int l = a.length ();
+  if (l != b.length ())
+    {
+      (*current_liboctave_error_handler)
+	("nonconformant array product attempted");
+      return MArray<T> ();
+    }
+
+  if (l == 0)
+    return MArray<T> ();
+
+  DO_VV_OP (*);
+  return MArray<T> (result, l);
+}
+
+template <class T>
+MArray<T>
+quotient (const MArray<T>& a, const MArray<T>& b)
+{
+  int l = a.length ();
+  if (l != b.length ())
+    {
+      (*current_liboctave_error_handler)
+	("nonconformant array quotient attempted");
+      return MArray<T> ();
+    }
+
+  if (l == 0)
+    return MArray<T> ();
+
+  DO_VV_OP (/);
+  return MArray<T> (result, l);
+}
+
+// Unary MArray ops.
+
+template <class T>
+MArray<T>
+operator - (const MArray<T>& a)
+{
+  NEG_V;
+  return MArray<T> (result, l);
+}
+
+/*
+ * Two dimensional array with math ops.
+ */
+
+template <class T>
+MArray2<T>::MArray2 (const MDiagArray<T>& a)
+  : Array2<T> (a.rows (), a.cols (), T (0))
+{
+  for (int i = 0; i < a.length (); i++)
+    elem (i, i) = a.elem (i, i);
+}
+
+// Element by element MArray2 by scalar ops.
+
+template <class T>
+MArray2<T>
+operator + (const MArray2<T>& a, const T& s)
+{
+  DO_VS_OP (+);
+  return MArray2<T> (result, a.rows (), a.cols ());
+}
+
+template <class T>
+MArray2<T>
+operator - (const MArray2<T>& a, const T& s)
+{
+  DO_VS_OP (-);
+  return MArray2<T> (result, a.rows (), a.cols ());
+}
+
+template <class T>
+MArray2<T>
+operator * (const MArray2<T>& a, const T& s)
+{
+  DO_VS_OP (*);
+  return MArray2<T> (result, a.rows (), a.cols ());
+}
+
+template <class T>
+MArray2<T>
+operator / (const MArray2<T>& a, const T& s)
+{
+  DO_VS_OP (/);
+  return MArray2<T> (result, a.rows (), a.cols ());
+}
+
+// Element by element scalar by MArray2 ops.
+
+template <class T>
+MArray2<T>
+operator + (const T& s, const MArray2<T>& a)
+{
+  DO_SV_OP (+);
+  return MArray2<T> (result, a.rows (), a.cols ());
+}
+
+template <class T>
+MArray2<T>
+operator - (const T& s, const MArray2<T>& a)
+{
+  DO_SV_OP (-);
+  return MArray2<T> (result, a.rows (), a.cols ());
+}
+
+template <class T>
+MArray2<T>
+operator * (const T& s, const MArray2<T>& a)
+{
+  DO_SV_OP (*);
+  return MArray2<T> (result, a.rows (), a.cols ());
+}
+
+template <class T>
+MArray2<T>
+operator / (const T& s, const MArray2<T>& a)
+{
+  DO_SV_OP (/);
+  return MArray2<T> (result, a.rows (), a.cols ());
+}
+
+// Element by element MArray2 by MArray2 ops.
+
+template <class T>
+MArray2<T>
+operator + (const MArray2<T>& a, const MArray2<T>& b)
+{
+  int r = a.rows ();
+  int c = a.cols ();
+  if (r != b.rows () || c != b.cols ())
+    {
+      (*current_liboctave_error_handler)
+	("nonconformant array addition attempted");
+      return MArray2<T> ();
+    }
+
+  if (r == 0 || c == 0)
+    return MArray2<T> ();
+
+  int l = a.length ();
+  DO_VV_OP (+);
+  return MArray2<T> (result, r, c);
+}
+
+template <class T>
+MArray2<T>
+operator - (const MArray2<T>& a, const MArray2<T>& b)
+{
+  int r = a.rows ();
+  int c = a.cols ();
+  if (r != b.rows () || c != b.cols ())
+    {
+      (*current_liboctave_error_handler)
+	("nonconformant array subtraction attempted");
+      return MArray2<T> ();
+    }
+
+  if (r == 0 || c == 0)
+    return MArray2<T> ();
+
+  int l = a.length ();
+  DO_VV_OP (-);
+  return MArray2<T> (result, r, c);
+}
+
+template <class T>
+MArray2<T>
+product (const MArray2<T>& a, const MArray2<T>& b)
+{
+  int r = a.rows ();
+  int c = a.cols ();
+  if (r != b.rows () || c != b.cols ())
+    {
+      (*current_liboctave_error_handler)
+	("nonconformant array product attempted");
+      return MArray2<T> ();
+    }
+
+  if (r == 0 || c == 0)
+    return MArray2<T> ();
+
+  int l = a.length ();
+  DO_VV_OP (*);
+  return MArray2<T> (result, r, c);
+}
+
+template <class T>
+MArray2<T>
+quotient (const MArray2<T>& a, const MArray2<T>& b)
+{
+  int r = a.rows ();
+  int c = a.cols ();
+  if (r != b.rows () || c != b.cols ())
+    {
+      (*current_liboctave_error_handler)
+	("nonconformant array quotient attempted");
+      return MArray2<T> ();
+    }
+
+  if (r == 0 || c == 0)
+    return MArray2<T> ();
+
+  int l = a.length ();
+  DO_VV_OP (/);
+  return MArray2<T> (result, r, c);
+}
+
+// Unary MArray2 ops.
+
+template <class T>
+MArray2<T>
+operator - (const MArray2<T>& a)
+{
+  NEG_V;
+  return MArray2<T> (result, a.rows (), a.cols ());
+}
+
+/*
+ * Two dimensional diagonal array with math ops.
+ */
+
+// Element by element MDiagArray by scalar ops.
+
+template <class T>
+MDiagArray<T>
+operator * (const MDiagArray<T>& a, const T& s)
+{
+  DO_VS_OP (*);
+  return MDiagArray<T> (result, a.rows (), a.cols ());
+}
+
+template <class T>
+MDiagArray<T>
+operator / (const MDiagArray<T>& a, const T& s)
+{
+  DO_VS_OP (/);
+  return MDiagArray<T> (result, a.rows (), a.cols ());
+}
+
+// Element by element scalar by MDiagArray ops.
+
+template <class T>
+MDiagArray<T>
+operator * (const T& s, const MDiagArray<T>& a)
+{
+  DO_SV_OP (*);
+  return MDiagArray<T> (result, a.rows (), a.cols ());
+}
+
+// Element by element MDiagArray by MDiagArray ops.
+
+template <class T>
+MDiagArray<T>
+operator + (const MDiagArray<T>& a, const MDiagArray<T>& b)
+{
+  int r = a.rows ();
+  int c = a.cols ();
+  if (r != b.rows () || c != b.cols ())
+    {
+      (*current_liboctave_error_handler)
+	("nonconformant diagonal array addition attempted");
+      return MDiagArray<T> ();
+    }
+
+  if (c == 0 || r == 0)
+    return MDiagArray<T> ();
+
+  int l = a.length ();
+  DO_VV_OP (+);
+  return MDiagArray<T> (result, r, c);
+}
+
+template <class T>
+MDiagArray<T>
+operator - (const MDiagArray<T>& a, const MDiagArray<T>& b)
+{
+  int r = a.rows ();
+  int c = a.cols ();
+  if (r != b.rows () || c != b.cols ())
+    {
+      (*current_liboctave_error_handler)
+	("nonconformant diagonal array subtraction attempted");
+      return MDiagArray<T> ();
+    }
+
+  if (c == 0 || r == 0)
+    return MDiagArray<T> ();
+
+  int l = a.length ();
+  DO_VV_OP (-);
+  return MDiagArray<T> (result, r, c);
+}
+
+template <class T>
+MDiagArray<T>
+product (const MDiagArray<T>& a, const MDiagArray<T>& b)
+{
+  int r = a.rows ();
+  int c = a.cols ();
+  if (r != b.rows () || c != b.cols ())
+    {
+      (*current_liboctave_error_handler)
+	("nonconformant diagonal array product attempted");
+      return MDiagArray<T> ();
+    }
+
+  if (c == 0 || r == 0)
+    return MDiagArray<T> ();
+
+  int l = a.length ();
+  DO_VV_OP (*);
+  return MDiagArray<T> (result, r, c);
+}
+
+// Unary MDiagArray ops.
+
+template <class T>
+MDiagArray<T>
+operator - (const MDiagArray<T>& a)
+{
+  NEG_V;
+  return MDiagArray<T> (result, a.rows (), a.cols ());
+}
+
+#undef DO_SV_OP
+#undef DO_VS_OP
+#undef DO_VV_OP
+#undef NEG_V
+
+#if 0
+#ifdef OCTAVE
+typedefMArray<double>      octave_mad_template_type;
+typedefMArray2<double>     octave_ma2d_template_type;
+typedefMDiagArray<double>  octave_mdad_template_type;
+
+#include <Complex.h>
+typedefMArray<Complex>     octave_mac_template_type;
+typedefMArray2<Complex>    octave_ma2c_template_type;
+typedefMDiagArray<Complex> octave_mdac_template_type;
+#endif
+#endif
+
+/*
+;;; Local Variables: ***
+;;; mode: C++ ***
+;;; page-delimiter: "^/\\*" ***
+;;; End: ***
+*/
new file mode 100644
--- /dev/null
+++ b/liboctave/MArray.h
@@ -0,0 +1,213 @@
+// Template array classes with like-type math ops          -*- C++ -*-
+/*
+
+Copyright (C) 1993 John W. Eaton
+
+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 2, 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, write to the Free
+Software Foundation, 675 Mass Ave, Cambridge, MA 02139, USA.
+
+*/
+
+#if !defined (_MArray_h)
+#define _MArray_h 1
+
+#include "Array.h"
+
+// Classes we declare.
+
+template <class T> class MArray;
+template <class T> class MArray2;
+template <class T> class MDiagArray;
+
+/*
+ * One dimensional array with math ops.
+ */
+
+template <class T>
+class MArray : public Array<T>
+{
+protected:
+
+  MArray (T *d, int l) : Array<T> (d, l) { }
+
+public:
+  
+  MArray (void) : Array<T> () { }
+  MArray (int n) : Array<T> (n) { }
+  MArray (int n, const T& val) : Array<T> (n, val) { }
+  MArray (const Array<T>& a) : Array<T> (a) { }
+  MArray (const MArray<T>& a) : Array<T> (a) { }
+
+  MArray<T>& operator = (const MArray<T>& a)
+    { return Array<T>::operator = (a); }
+
+// Element by element MArray by scalar ops.
+
+  MArray<T>& operator += (const T& s) { /* XXX FIXME XXX */ }
+  MArray<T>& operator -= (const T& s) { /* XXX FIXME XXX */ }
+
+// Element by element MArray by MArray ops.
+
+  MArray<T>& operator += (const MArray<T>& a) { /* XXX FIXME XXX */ }
+  MArray<T>& operator -= (const MArray<T>& a) { /* XXX FIXME XXX */ }
+
+// Element by element MArray by scalar ops.
+
+  friend MArray<T> operator + (const MArray<T>& a, const T& s);
+  friend MArray<T> operator - (const MArray<T>& a, const T& s);
+  friend MArray<T> operator * (const MArray<T>& a, const T& s);
+  friend MArray<T> operator / (const MArray<T>& a, const T& s);
+
+// Element by element scalar by MArray ops.
+
+  friend MArray<T> operator + (const T& s, const MArray<T>& a);
+  friend MArray<T> operator - (const T& s, const MArray<T>& a);
+  friend MArray<T> operator * (const T& s, const MArray<T>& a);
+  friend MArray<T> operator / (const T& s, const MArray<T>& a);
+
+// Element by element MArray by MArray ops.
+
+  friend MArray<T> operator + (const MArray<T>& a, const MArray<T>& b);
+
+  friend MArray<T> operator - (const MArray<T>& a, const MArray<T>& b);
+
+  friend MArray<T> product (const MArray<T>& a, const MArray<T>& b);
+  friend MArray<T> quotient (const MArray<T>& a, const MArray<T>& b);
+
+  friend MArray<T> operator - (const MArray<T>& a);
+};
+
+/*
+ * Two dimensional array with math ops.
+ */
+
+template <class T>
+class MArray2 : public Array2<T>
+{
+protected:
+
+  MArray2 (T *d, int n, int m) : Array2<T> (d, n, m) { }
+
+public:
+
+  MArray2 (void) : Array2<T> () { }
+  MArray2 (int n, int m) : Array2<T> (n, m) { }
+  MArray2 (int n, int m, const T& val) : Array2<T> (n, m, val) { }
+  MArray2 (const Array2<T>& a) : Array2<T> (a) { }
+  MArray2 (const MArray2<T>& a) : Array2<T> (a) { }
+  MArray2 (const MDiagArray<T>& a);
+
+  MArray2<T>& operator = (const MArray2<T>& a)
+    { return Array2<T>::operator = (a); }
+
+// Element by element MArray2 by scalar ops.
+
+  MArray2<T>& operator += (const T& s) { /* XXX FIXME XXX */ }
+  MArray2<T>& operator -= (const T& s) { /* XXX FIXME XXX */ }
+
+// Element by element MArray2 by MArray2 ops.
+
+  MArray2<T>& operator += (const MArray2<T>& a) { /* XXX FIXME XXX */ }
+  MArray2<T>& operator -= (const MArray2<T>& a) { /* XXX FIXME XXX */ }
+
+// Element by element MArray2 by scalar ops.
+
+  friend MArray2<T> operator + (const MArray2<T>& a, const T& s);
+  friend MArray2<T> operator - (const MArray2<T>& a, const T& s);
+  friend MArray2<T> operator * (const MArray2<T>& a, const T& s);
+  friend MArray2<T> operator / (const MArray2<T>& a, const T& s);
+
+// Element by element scalar by MArray2 ops.
+
+  friend MArray2<T> operator + (const T& s, const MArray2<T>& a);
+  friend MArray2<T> operator - (const T& s, const MArray2<T>& a);
+  friend MArray2<T> operator * (const T& s, const MArray2<T>& a);
+  friend MArray2<T> operator / (const T& s, const MArray2<T>& a);
+
+// Element by element MArray2 by MArray2 ops.
+
+  friend MArray2<T> operator + (const MArray2<T>& a, const MArray2<T>& b);
+  friend MArray2<T> operator - (const MArray2<T>& a, const MArray2<T>& b);
+
+  friend MArray2<T> product (const MArray2<T>& a, const MArray2<T>& b);
+  friend MArray2<T> quotient (const MArray2<T>& a, const MArray2<T>& b);
+
+  friend MArray2<T> operator - (const MArray2<T>& a);
+};
+
+/*
+ * Two dimensional diagonal array with math ops.
+ */
+
+template <class T>
+class MDiagArray : public DiagArray<T>
+{
+protected:
+
+  MDiagArray (T *d, int r, int c) : DiagArray<T> (d, r, c) { }
+
+public:
+  
+  MDiagArray (void) : DiagArray<T> () { }
+  MDiagArray (int n) : DiagArray<T> (n) { }
+  MDiagArray (int n, const T& val) : DiagArray<T> (n, val) { }
+  MDiagArray (int r, int c) : DiagArray<T> (r, c) { }
+  MDiagArray (int r, int c, const T& val) : DiagArray<T> (r, c, val) { }
+  MDiagArray (const DiagArray<T>& a) : DiagArray<T> (a) { }
+  MDiagArray (const MDiagArray<T>& a) : DiagArray<T> (a) { }
+  MDiagArray (const MArray<T>& a) : DiagArray<T> (a) { }
+
+  MDiagArray<T>& operator = (const MDiagArray<T>& a)
+    { return DiagArray<T>::operator = (a); }
+
+  MDiagArray<T>& operator += (const MDiagArray<T>& a) { /* XXX FIXME XXX */ }
+  MDiagArray<T>& operator -= (const MDiagArray<T>& a) { /* XXX FIXME XXX */ }
+
+// Element by element MDiagArray by scalar ops.
+
+  friend MDiagArray<T> operator * (const MDiagArray<T>& a, const T& s);
+  friend MDiagArray<T> operator / (const MDiagArray<T>& a, const T& s);
+
+// Element by element scalar by MDiagArray ops.
+
+  friend MDiagArray<T> operator * (const T& s, const MDiagArray<T>& a);
+
+// Element by element MDiagArray by MDiagArray ops.
+
+  friend MDiagArray<T> operator + (const MDiagArray<T>& a,
+				   const MDiagArray<T>& b); 
+
+  friend MDiagArray<T> operator - (const MDiagArray<T>& a,
+				   const MDiagArray<T>& b);
+
+  friend MDiagArray<T> product (const MDiagArray<T>& a,
+				const MDiagArray<T>& b);
+
+  friend MDiagArray<T> operator - (const MDiagArray<T>& a);
+};
+
+#ifdef __GNUG__
+#include "MArray.cc"
+#endif
+
+#endif
+
+/*
+;;; Local Variables: ***
+;;; mode: C++ ***
+;;; page-delimiter: "^/\\*" ***
+;;; End: ***
+*/
new file mode 100644
--- /dev/null
+++ b/liboctave/mx-kludge.cc
@@ -0,0 +1,501 @@
+// kludge.cc                                             -*- C++ -*-
+/*
+
+Copyright (C) 1992, 1993 John W. Eaton
+
+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 2, 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, write to the Free
+Software Foundation, 675 Mass Ave, Cambridge, MA 02139, USA.
+
+*/
+
+// Nothing like a little CPP abuse to brighten everyone's day.  Would
+// have been nice to do this with template functions but as of 2.5.x,
+// g++ seems to fail in various ways, either not resolving general
+// template functions, or not instatiating non-member template
+// functions.
+//
+// When templates work more reliably in g++, this will be replaced by
+// the MArray class.
+
+#define DO_VS_OP(OP) \
+  int l = a.length (); \
+  TYPE *result = 0; \
+  if (l > 0) \
+    { \
+      result = new TYPE [l]; \
+      const TYPE *x = a.data (); \
+      for (int i = 0; i < l; i++) \
+	result[i] = x[i] OP s; \
+    }
+
+#define DO_SV_OP(OP) \
+  int l = a.length (); \
+  TYPE *result = 0; \
+  if (l > 0) \
+    { \
+      result = new TYPE [l]; \
+      const TYPE *x = a.data (); \
+      for (int i = 0; i < l; i++) \
+	result[i] = s OP x[i]; \
+    }
+
+#define DO_VV_OP(OP) \
+  TYPE *result = 0; \
+  if (l > 0) \
+    { \
+      result = new TYPE [l]; \
+      const TYPE *x = a.data (); \
+      const TYPE *y = b.data (); \
+      for (int i = 0; i < l; i++) \
+	result[i] = x[i] OP y[i]; \
+    }
+
+#define NEG_V \
+  int l = a.length (); \
+  TYPE *result = 0; \
+  if (l > 0) \
+    { \
+      result = new TYPE [l]; \
+      const TYPE *x = a.data (); \
+      for (int i = 0; i < l; i++) \
+	result[i] = -x[i]; \
+    }
+
+#ifdef KLUDGE_VECTORS
+
+/*
+ * Like type operations for vectors.
+ */
+
+// Element by element vector by scalar ops.
+
+KL_VEC_TYPE
+operator + (const KL_VEC_TYPE& a, const TYPE& s)
+{
+  DO_VS_OP (+);
+  return KL_VEC_TYPE (result, l);
+}
+
+KL_VEC_TYPE
+operator - (const KL_VEC_TYPE& a, const TYPE& s)
+{
+  DO_VS_OP (-);
+  return KL_VEC_TYPE (result, l);
+}
+
+KL_VEC_TYPE
+operator * (const KL_VEC_TYPE& a, const TYPE& s)
+{
+  DO_VS_OP (*);
+  return KL_VEC_TYPE (result, l);
+}
+
+KL_VEC_TYPE
+operator / (const KL_VEC_TYPE& a, const TYPE& s)
+{
+  DO_VS_OP (/);
+  return KL_VEC_TYPE (result, l);
+}
+
+// Element by element scalar by vector ops.
+
+KL_VEC_TYPE
+operator + (const TYPE& s, const KL_VEC_TYPE& a)
+{
+  DO_SV_OP (+);
+  return KL_VEC_TYPE (result, l);
+}
+
+KL_VEC_TYPE
+operator - (const TYPE& s, const KL_VEC_TYPE& a)
+{
+  DO_SV_OP (-);
+  return KL_VEC_TYPE (result, l);
+}
+
+KL_VEC_TYPE
+operator * (const TYPE& s, const KL_VEC_TYPE& a)
+{
+  DO_SV_OP (*);
+  return KL_VEC_TYPE (result, l);
+}
+
+KL_VEC_TYPE
+operator / (const TYPE& s, const KL_VEC_TYPE& a)
+{
+  DO_SV_OP (/);
+  return KL_VEC_TYPE (result, l);
+}
+
+// Element by element vector by vector ops.
+
+KL_VEC_TYPE
+operator + (const KL_VEC_TYPE& a, const KL_VEC_TYPE& b)
+{
+  int l = a.length ();
+  if (l != b.length ())
+    {
+      (*current_liboctave_error_handler)
+	("nonconformant array addition attempted");
+      return KL_VEC_TYPE ();
+    }
+
+  if (l == 0)
+    return KL_VEC_TYPE ();
+
+  DO_VV_OP (+);
+  return KL_VEC_TYPE (result, l);
+}
+
+KL_VEC_TYPE
+operator - (const KL_VEC_TYPE& a, const KL_VEC_TYPE& b)
+{
+  int l = a.length ();
+  if (l != b.length ())
+    {
+      (*current_liboctave_error_handler)
+	("nonconformant array subtraction attempted");
+      return KL_VEC_TYPE ();
+    }
+
+  if (l == 0)
+    return KL_VEC_TYPE ();
+
+  DO_VV_OP (-);
+  return KL_VEC_TYPE (result, l);
+}
+
+KL_VEC_TYPE
+product (const KL_VEC_TYPE& a, const KL_VEC_TYPE& b)
+{
+  int l = a.length ();
+  if (l != b.length ())
+    {
+      (*current_liboctave_error_handler)
+	("nonconformant array product attempted");
+      return KL_VEC_TYPE ();
+    }
+
+  if (l == 0)
+    return KL_VEC_TYPE ();
+
+  DO_VV_OP (*);
+  return KL_VEC_TYPE (result, l);
+}
+
+KL_VEC_TYPE
+quotient (const KL_VEC_TYPE& a, const KL_VEC_TYPE& b)
+{
+  int l = a.length ();
+  if (l != b.length ())
+    {
+      (*current_liboctave_error_handler)
+	("nonconformant array quotient attempted");
+      return KL_VEC_TYPE ();
+    }
+
+  if (l == 0)
+    return KL_VEC_TYPE ();
+
+  DO_VV_OP (/);
+  return KL_VEC_TYPE (result, l);
+}
+
+// Unary MArray ops.
+
+KL_VEC_TYPE
+operator - (const KL_VEC_TYPE& a)
+{
+  NEG_V;
+  return KL_VEC_TYPE (result, l);
+}
+
+#endif
+
+#ifdef KLUDGE_MATRICES
+
+/*
+ * Like type operations for matrices
+ */
+
+// Element by element matrix by scalar ops.
+
+KL_MAT_TYPE
+operator + (const KL_MAT_TYPE& a, const TYPE& s)
+{
+  DO_VS_OP (+);
+  return KL_MAT_TYPE (result, a.rows (), a.cols ());
+}
+
+KL_MAT_TYPE
+operator - (const KL_MAT_TYPE& a, const TYPE& s)
+{
+  DO_VS_OP (-);
+  return KL_MAT_TYPE (result, a.rows (), a.cols ());
+}
+
+KL_MAT_TYPE
+operator * (const KL_MAT_TYPE& a, const TYPE& s)
+{
+  DO_VS_OP (*);
+  return KL_MAT_TYPE (result, a.rows (), a.cols ());
+}
+
+KL_MAT_TYPE
+operator / (const KL_MAT_TYPE& a, const TYPE& s)
+{
+  DO_VS_OP (/);
+  return KL_MAT_TYPE (result, a.rows (), a.cols ());
+}
+
+// Element by element scalar by matrix ops.
+
+KL_MAT_TYPE
+operator + (const TYPE& s, const KL_MAT_TYPE& a)
+{
+  DO_SV_OP (+);
+  return KL_MAT_TYPE (result, a.rows (), a.cols ());
+}
+
+KL_MAT_TYPE
+operator - (const TYPE& s, const KL_MAT_TYPE& a)
+{
+  DO_SV_OP (-);
+  return KL_MAT_TYPE (result, a.rows (), a.cols ());
+}
+
+KL_MAT_TYPE
+operator * (const TYPE& s, const KL_MAT_TYPE& a)
+{
+  DO_SV_OP (*);
+  return KL_MAT_TYPE (result, a.rows (), a.cols ());
+}
+
+KL_MAT_TYPE
+operator / (const TYPE& s, const KL_MAT_TYPE& a)
+{
+  DO_SV_OP (/);
+  return KL_MAT_TYPE (result, a.rows (), a.cols ());
+}
+
+// Element by element matrix by matrix ops.
+
+KL_MAT_TYPE
+operator + (const KL_MAT_TYPE& a, const KL_MAT_TYPE& b)
+{
+  int r = a.rows ();
+  int c = a.cols ();
+  if (r != b.rows () || c != b.cols ())
+    {
+      (*current_liboctave_error_handler)
+	("nonconformant array addition attempted");
+      return KL_MAT_TYPE ();
+    }
+
+  if (r == 0 || c == 0)
+    return KL_MAT_TYPE ();
+
+  int l = a.length ();
+  DO_VV_OP (+);
+  return KL_MAT_TYPE (result, r, c);
+}
+
+KL_MAT_TYPE
+operator - (const KL_MAT_TYPE& a, const KL_MAT_TYPE& b)
+{
+  int r = a.rows ();
+  int c = a.cols ();
+  if (r != b.rows () || c != b.cols ())
+    {
+      (*current_liboctave_error_handler)
+	("nonconformant array subtraction attempted");
+      return KL_MAT_TYPE ();
+    }
+
+  if (r == 0 || c == 0)
+    return KL_MAT_TYPE ();
+
+  int l = a.length ();
+  DO_VV_OP (-);
+  return KL_MAT_TYPE (result, r, c);
+}
+
+KL_MAT_TYPE
+product (const KL_MAT_TYPE& a, const KL_MAT_TYPE& b)
+{
+  int r = a.rows ();
+  int c = a.cols ();
+  if (r != b.rows () || c != b.cols ())
+    {
+      (*current_liboctave_error_handler)
+	("nonconformant array product attempted");
+      return KL_MAT_TYPE ();
+    }
+
+  if (r == 0 || c == 0)
+    return KL_MAT_TYPE ();
+
+  int l = a.length ();
+  DO_VV_OP (*);
+  return KL_MAT_TYPE (result, r, c);
+}
+
+KL_MAT_TYPE
+quotient (const KL_MAT_TYPE& a, const KL_MAT_TYPE& b)
+{
+  int r = a.rows ();
+  int c = a.cols ();
+  if (r != b.rows () || c != b.cols ())
+    {
+      (*current_liboctave_error_handler)
+	("nonconformant array quotient attempted");
+      return KL_MAT_TYPE ();
+    }
+
+  if (r == 0 || c == 0)
+    return KL_MAT_TYPE ();
+
+  int l = a.length ();
+  DO_VV_OP (/);
+  return KL_MAT_TYPE (result, r, c);
+}
+
+// Unary matrix ops.
+
+KL_MAT_TYPE
+operator - (const KL_MAT_TYPE& a)
+{
+  NEG_V;
+  return KL_MAT_TYPE (result, a.rows (), a.cols ());
+}
+
+#endif
+
+#ifdef KLUDGE_DIAG_MATRICES
+
+/*
+ * Like type operations for diagonal matrices.
+ */
+
+// Element by element MDiagArray by scalar ops.
+
+KL_DMAT_TYPE
+operator * (const KL_DMAT_TYPE& a, const TYPE& s)
+{
+  DO_VS_OP (*);
+  return KL_DMAT_TYPE (result, a.rows (), a.cols ());
+}
+
+KL_DMAT_TYPE
+operator / (const KL_DMAT_TYPE& a, const TYPE& s)
+{
+  DO_VS_OP (/);
+  return KL_DMAT_TYPE (result, a.rows (), a.cols ());
+}
+
+// Element by element scalar by MDiagArray ops.
+
+KL_DMAT_TYPE
+operator * (const TYPE& s, const KL_DMAT_TYPE& a)
+{
+  DO_SV_OP (*);
+  return KL_DMAT_TYPE (result, a.rows (), a.cols ());
+}
+
+// Element by element MDiagArray by MDiagArray ops.
+
+KL_DMAT_TYPE
+operator + (const KL_DMAT_TYPE& a, const KL_DMAT_TYPE& b)
+{
+  int r = a.rows ();
+  int c = a.cols ();
+  if (r != b.rows () || c != b.cols ())
+    {
+      (*current_liboctave_error_handler)
+	("nonconformant diagonal array addition attempted");
+      return KL_DMAT_TYPE ();
+    }
+
+  if (c == 0 || r == 0)
+    return KL_DMAT_TYPE ();
+
+  int l = a.length ();
+  DO_VV_OP (+);
+  return KL_DMAT_TYPE (result, r, c);
+}
+
+KL_DMAT_TYPE
+operator - (const KL_DMAT_TYPE& a, const KL_DMAT_TYPE& b)
+{
+  int r = a.rows ();
+  int c = a.cols ();
+  if (r != b.rows () || c != b.cols ())
+    {
+      (*current_liboctave_error_handler)
+	("nonconformant diagonal array subtraction attempted");
+      return KL_DMAT_TYPE ();
+    }
+
+  if (c == 0 || r == 0)
+    return KL_DMAT_TYPE ();
+
+  int l = a.length ();
+  DO_VV_OP (-);
+  return KL_DMAT_TYPE (result, r, c);
+}
+
+KL_DMAT_TYPE
+product (const KL_DMAT_TYPE& a, const KL_DMAT_TYPE& b)
+{
+  int r = a.rows ();
+  int c = a.cols ();
+  if (r != b.rows () || c != b.cols ())
+    {
+      (*current_liboctave_error_handler)
+	("nonconformant diagonal array product attempted");
+      return KL_DMAT_TYPE ();
+    }
+
+  if (c == 0 || r == 0)
+    return KL_DMAT_TYPE ();
+
+  int l = a.length ();
+  DO_VV_OP (*);
+  return KL_DMAT_TYPE (result, r, c);
+}
+
+// Unary MDiagArray ops.
+
+KL_DMAT_TYPE
+operator - (const KL_DMAT_TYPE& a)
+{
+  NEG_V;
+  return KL_DMAT_TYPE (result, a.rows (), a.cols ());
+}
+
+#endif
+
+#undef DO_SV_OP
+#undef DO_VS_OP
+#undef DO_VV_OP
+#undef NEG_V
+
+/*
+;;; Local Variables: ***
+;;; mode: C++ ***
+;;; page-delimiter: "^/\\*" ***
+;;; End: ***
+*/
new file mode 100644
--- /dev/null
+++ b/liboctave/mx-kludge.h
@@ -0,0 +1,139 @@
+// kludge.h                                             -*- C++ -*-
+/*
+
+Copyright (C) 1992, 1993 John W. Eaton
+
+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 2, 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, write to the Free
+Software Foundation, 675 Mass Ave, Cambridge, MA 02139, USA.
+
+*/
+
+// Nothing like a little CPP abuse to brighten everyone's day.  Would
+// have been nice to do this with template functions but as of 2.5.x,
+// g++ seems to fail in various ways, either not resolving general
+// template functions, or not instatiating non-member template
+// functions.
+//
+// When templates work more reliably in g++, this will be replaced by
+// the MArray class.
+
+#ifdef KLUDGE_VECTORS
+
+/*
+ * Like type operations for vectors.
+ */
+
+// Element by element vector by scalar ops.
+
+friend KL_VEC_TYPE operator + (const KL_VEC_TYPE& a, const TYPE& s);
+friend KL_VEC_TYPE operator - (const KL_VEC_TYPE& a, const TYPE& s);
+friend KL_VEC_TYPE operator * (const KL_VEC_TYPE& a, const TYPE& s);
+friend KL_VEC_TYPE operator / (const KL_VEC_TYPE& a, const TYPE& s);
+
+// Element by element scalar by vector ops.
+
+friend KL_VEC_TYPE operator + (const TYPE& s, const KL_VEC_TYPE& a);
+friend KL_VEC_TYPE operator - (const TYPE& s, const KL_VEC_TYPE& a);
+friend KL_VEC_TYPE operator * (const TYPE& s, const KL_VEC_TYPE& a);
+friend KL_VEC_TYPE operator / (const TYPE& s, const KL_VEC_TYPE& a);
+
+// Element by element vector by vector ops.
+
+friend KL_VEC_TYPE operator + (const KL_VEC_TYPE& a, const KL_VEC_TYPE& b);
+friend KL_VEC_TYPE operator - (const KL_VEC_TYPE& a, const KL_VEC_TYPE& b);
+
+friend KL_VEC_TYPE product (const KL_VEC_TYPE& a, const KL_VEC_TYPE& b);
+friend KL_VEC_TYPE quotient (const KL_VEC_TYPE& a, const KL_VEC_TYPE& b);
+
+// Unary MArray ops.
+
+friend KL_VEC_TYPE operator - (const KL_VEC_TYPE& a);
+
+#endif
+
+#ifdef KLUDGE_MATRICES
+
+/*
+ * Like type operations for matrices
+ */
+
+// Element by element matrix by scalar ops.
+
+friend KL_MAT_TYPE operator + (const KL_MAT_TYPE& a, const TYPE& s);
+friend KL_MAT_TYPE operator - (const KL_MAT_TYPE& a, const TYPE& s);
+friend KL_MAT_TYPE operator * (const KL_MAT_TYPE& a, const TYPE& s);
+friend KL_MAT_TYPE operator / (const KL_MAT_TYPE& a, const TYPE& s);
+
+// Element by element scalar by matrix ops.
+
+friend KL_MAT_TYPE operator + (const TYPE& s, const KL_MAT_TYPE& a);
+friend KL_MAT_TYPE operator - (const TYPE& s, const KL_MAT_TYPE& a);
+friend KL_MAT_TYPE operator * (const TYPE& s, const KL_MAT_TYPE& a);
+friend KL_MAT_TYPE operator / (const TYPE& s, const KL_MAT_TYPE& a);
+
+// Element by element matrix by matrix ops.
+
+friend KL_MAT_TYPE operator + (const KL_MAT_TYPE& a, const KL_MAT_TYPE& b);
+friend KL_MAT_TYPE operator - (const KL_MAT_TYPE& a, const KL_MAT_TYPE& b);
+
+friend KL_MAT_TYPE product (const KL_MAT_TYPE& a, const KL_MAT_TYPE& b);
+friend KL_MAT_TYPE quotient (const KL_MAT_TYPE& a, const KL_MAT_TYPE& b);
+
+// Unary matrix ops.
+
+friend KL_MAT_TYPE operator - (const KL_MAT_TYPE& a);
+
+#endif
+
+#ifdef KLUDGE_DIAG_MATRICES
+
+/*
+ * Like type operations for diagonal matrices.
+ */
+
+// Element by element MDiagArray by scalar ops.
+
+friend KL_DMAT_TYPE operator * (const KL_DMAT_TYPE& a, const TYPE& s);
+friend KL_DMAT_TYPE operator / (const KL_DMAT_TYPE& a, const TYPE& s);
+
+// Element by element scalar by MDiagArray ops.
+
+friend KL_DMAT_TYPE operator * (const TYPE& s, const KL_DMAT_TYPE& a);
+
+// Element by element MDiagArray by MDiagArray ops.
+
+friend KL_DMAT_TYPE operator + (const KL_DMAT_TYPE& a, const KL_DMAT_TYPE& b);
+friend KL_DMAT_TYPE operator - (const KL_DMAT_TYPE& a, const KL_DMAT_TYPE& b);
+
+friend KL_DMAT_TYPE product (const KL_DMAT_TYPE& a, const KL_DMAT_TYPE& b);
+
+// Unary MDiagArray ops.
+
+friend KL_DMAT_TYPE operator - (const KL_DMAT_TYPE& a);
+
+#endif
+
+#undef DO_SV_OP
+#undef DO_VS_OP
+#undef DO_VV_OP
+#undef NEG_V
+
+/*
+;;; Local Variables: ***
+;;; mode: C++ ***
+;;; page-delimiter: "^/\\*" ***
+;;; End: ***
+*/