new file mode 100644
--- /dev/null
+++ b/src/f-eval.cc
@@ -0,0 +1,124 @@
+// f-eval.cc                                           -*- C++ -*-
+/*
+
+Copyright (C) 1994 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 "unwind-prot.h"
+#include "tree-const.h"
+#include "variables.h"
+#include "octave.h"
+#include "input.h"
+#include "parse.h"
+#include "lex.h"
+#include "f-eval.h"
+
+Octave_object
+feval (const Octave_object& args, int nargout)
+{
+// Assumes that we have been given the correct number of arguments.
+
+  Octave_object retval;
+
+  tree_fvc *fcn = is_valid_function (args(1), "feval", 1);
+  if (fcn != (tree_fvc *) NULL)
+    {
+      int nargin = args.length () - 1;
+      Octave_object tmp_args (nargin);
+      for (int i = 0; i < nargin; i++)
+	tmp_args(i) = args(i+1);
+      retval = fcn->eval (0, nargout, tmp_args);
+    }
+
+  return retval;
+}
+
+tree_constant
+eval_string (const char *string, int print, int ans_assign,
+	     int& parse_status)
+{
+  begin_unwind_frame ("eval_string");
+
+  unwind_protect_int (get_input_from_eval_string);
+  unwind_protect_ptr (global_command);
+  unwind_protect_ptr (current_eval_string);
+
+  get_input_from_eval_string = 1;
+  current_eval_string = string;
+
+  YY_BUFFER_STATE old_buf = current_buffer ();
+  YY_BUFFER_STATE new_buf = create_buffer ((FILE *) NULL);
+
+  add_unwind_protect (restore_input_buffer, (void *) old_buf);
+  add_unwind_protect (delete_input_buffer, (void *) new_buf);
+
+  switch_to_buffer (new_buf);
+
+  unwind_protect_ptr (curr_sym_tab);
+
+  reset_parser ();
+
+  parse_status = yyparse ();
+
+// Important to reset the idea of where input is coming from before
+// trying to eval the command we just parsed -- it might contain the
+// name of an function file that still needs to be parsed!
+
+  tree *command = global_command;
+
+  run_unwind_frame ("eval_string");
+
+  tree_constant retval;
+
+  if (parse_status == 0 && command != NULL_TREE)
+    {
+      retval = command->eval (print);
+      delete command;
+    }
+
+  return retval;
+}
+
+tree_constant
+eval_string (const tree_constant& arg, int& parse_status)
+{
+  if (! arg.is_string_type ())
+    {
+      error ("eval: expecting string argument");
+      return -1;
+    }
+
+  char *string = arg.string_value ();
+
+// Yes Virginia, we always print here...
+
+  return eval_string (string, 1, 1, parse_status);
+}
+
+/*
+;;; Local Variables: ***
+;;; mode: C++ ***
+;;; page-delimiter: "^/\\*" ***
+;;; End: ***
+*/

new file mode 100644
--- /dev/null
+++ b/src/f-fill.cc
@@ -0,0 +1,175 @@
+// f-fill.cc                                           -*- C++ -*-
+/*
+
+Copyright (C) 1994 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 "tree-const.h"
+#include "user-prefs.h"
+#include "utils.h"
+#include "error.h"
+#include "f-fill.h"
+
+#ifndef MIN
+#define MIN(a,b) ((a) < (b) ? (a) : (b))
+#endif
+
+static void
+check_dimensions (int& nr, int& nc, const char *warnfor)
+{
+  if (nr < 0 || nc < 0)
+    {
+      if (user_pref.treat_neg_dim_as_zero)
+	nr = nc = 0;
+      else
+	error ("%s: can't create a matrix with negative dimensions",
+	       warnfor);
+    }
+}
+
+static void
+get_dimensions (const tree_constant& a, const char *warn_for,
+		int& nr, int& nc)
+{
+  tree_constant tmpa = a.make_numeric ();
+
+  if (tmpa.is_scalar_type ())
+    {
+      double tmp = tmpa.double_value ();
+      nr = nc = NINT (tmp);
+    }
+  else
+    {
+      nr = tmpa.rows ();
+      nc = tmpa.columns ();
+
+      if ((nr == 1 && nc == 2) || (nr == 2 && nc == 1))
+	{
+	  ColumnVector v = tmpa.to_vector ();
+
+	  nr = NINT (v.elem (0));
+	  nc = NINT (v.elem (1));
+	}
+      else
+	warning ("%s (A): use %s (size (A)) instead", warn_for, warn_for);
+    }
+
+  check_dimensions (nr, nc, warn_for); // May set error_state.
+}
+
+static void
+get_dimensions (const tree_constant& a, const tree_constant& b,
+		const char *warn_for, int& nr, int& nc)
+{
+  tree_constant tmpa = a.make_numeric ();
+  tree_constant tmpb = b.make_numeric ();
+
+  if (tmpa.is_scalar_type () && tmpb.is_scalar_type ())
+    {
+      nr = NINT (tmpa.double_value ());
+      nc = NINT (tmpb.double_value ());
+
+      check_dimensions (nr, nc, warn_for); // May set error_state.
+    }
+  else
+    error ("%s: expecting two scalar arguments", warn_for);
+}
+
+tree_constant
+fill_matrix (const tree_constant& a, double val, const char *warn_for)
+{
+  int nr, nc;
+  get_dimensions (a, warn_for, nr, nc);
+
+  if (error_state)
+    return  tree_constant ();
+
+  Matrix m (nr, nc, val);
+
+  return m;
+}
+
+tree_constant
+fill_matrix (const tree_constant& a, const tree_constant& b,
+	     double val, const char *warn_for)
+{
+  int nr, nc;
+  get_dimensions (a, b, warn_for, nr, nc); // May set error_state.
+
+  if (error_state)
+    return tree_constant ();
+
+  Matrix m (nr, nc, val);
+
+  return m;
+}
+
+tree_constant
+identity_matrix (const tree_constant& a)
+{
+  int nr, nc;
+  get_dimensions (a, "eye", nr, nc); // May set error_state.
+
+  if (error_state)
+    return tree_constant ();
+
+  Matrix m (nr, nc, 0.0);
+
+  if (nr > 0 && nc > 0)
+    {
+      int n = MIN (nr, nc);
+      for (int i = 0; i < n; i++)
+	m.elem (i, i) = 1.0;
+    }
+
+  return m;
+}
+
+tree_constant
+identity_matrix (const tree_constant& a, const tree_constant& b)
+{
+  int nr, nc;
+  get_dimensions (a, b, "eye", nr, nc);  // May set error_state.
+
+  if (error_state)
+    return tree_constant ();
+
+  Matrix m (nr, nc, 0.0);
+
+  if (nr > 0 && nc > 0)
+    {
+      int n = MIN (nr, nc);
+      for (int i = 0; i < n; i++)
+	m.elem (i, i) = 1.0;
+    }
+
+  return m;
+}
+
+/*
+;;; Local Variables: ***
+;;; mode: C++ ***
+;;; page-delimiter: "^/\\*" ***
+;;; End: ***
+*/

new file mode 100644
--- /dev/null
+++ b/src/f-input.cc
@@ -0,0 +1,135 @@
+// f-input.cc                                           -*- C++ -*-
+/*
+
+Copyright (C) 1994 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 <string.h>
+
+#include "tree-const.h"
+#include "octave-hist.h"
+#include "pager.h"
+#include "input.h"
+#include "f-eval.h"
+#include "f-input.h"
+
+static int
+match_sans_spaces (const char *standard, const char *test)
+{
+  const char *tp = test;
+  while (*tp == ' ' || *tp == '\t')
+    tp++;
+
+  const char *ep = test + strlen (test) - 1;
+  while (*ep == ' ' || *ep == '\t')
+    ep--;
+
+  int len = ep - tp + 1;
+
+  return (strncmp (standard, tp, len) == 0);
+}
+
+tree_constant
+get_user_input (const Octave_object& args, int nargout, int debug = 0)
+{
+  tree_constant retval;
+
+  int nargin = args.length ();
+
+  int read_as_string = 0;
+
+  if (nargin == 3)
+    {
+      if (args(2).is_string_type ()
+	  && strcmp ("s", args(2).string_value ()) == 0)
+	read_as_string++;
+      else
+	{
+	  error ("input: unrecognized second argument");
+	  return retval;
+	}
+    }
+
+  char *prompt = "debug> ";
+  if (nargin > 1)
+   {
+      if (args(1).is_string_type ())
+	prompt = args(1).string_value ();
+      else
+	{
+	  error ("input: unrecognized argument");
+	  return retval;
+	}
+    }
+
+ again:
+
+  flush_output_to_pager ();
+
+  char *input_buf = gnu_readline (prompt);
+
+  if (input_buf != (char *) NULL)
+    {
+      if (input_buf)
+	maybe_save_history (input_buf);
+
+      int len = strlen (input_buf);
+
+      if (len < 1)
+	{
+	  if (debug)
+	    goto again;
+	  else
+	    return retval;
+	}
+
+      if (match_sans_spaces ("exit", input_buf)
+	  || match_sans_spaces ("quit", input_buf)
+	  || match_sans_spaces ("return", input_buf))
+	return tree_constant ();
+      else if (read_as_string)
+	retval = input_buf;
+      else
+	{
+	  int parse_status = 0;
+	  retval = eval_string (input_buf, 0, 0, parse_status);
+	  if (debug && retval.is_defined ())
+	    retval.eval (1);
+	}
+    }
+  else
+    error ("input: reading user-input failed!");
+
+  if (debug)
+    goto again;
+
+  return retval;
+}
+
+/*
+;;; Local Variables: ***
+;;; mode: C++ ***
+;;; page-delimiter: "^/\\*" ***
+;;; End: ***
+*/

new file mode 100644
--- /dev/null
+++ b/src/find.cc
@@ -0,0 +1,212 @@
+// f-find.cc                                           -*- C++ -*-
+/*
+
+Copyright (C) 1994 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 "tree-const.h"
+#include "error.h"
+#include "f-find.h"
+
+static Octave_object
+find_to_fortran_idx (const ColumnVector i_idx, const ColumnVector j_idx,
+		     const tree_constant& val, int nr, int nc, int nargout)
+{
+  Octave_object retval (nargout);
+
+  switch (nargout)
+    {
+    case 1:
+      {
+	int count = i_idx.length ();
+	ColumnVector tmp (count);
+	for (int i = 0; i < count; i++)
+	  tmp (i) = nr * (j_idx (i) - 1.0) + i_idx (i);
+	retval(0) = tree_constant (tmp, 1);
+// If you want this to work more like Matlab, use the following line
+// instead of the previous one.
+//	retval(0) = tree_constant (tmp, (nr != 1));
+      }
+      break;
+    case 3:
+      retval(2) = val;
+    case 2:
+      retval(0) = tree_constant (i_idx, 1);
+// If you want this to work more like Matlab, use the following line
+// instead of the previous one.
+//    retval(0) = tree_constant (i_idx, (nr != 1));
+      retval(1) = tree_constant (j_idx, 1);
+      break;
+    default:
+      panic_impossible ();
+      break;
+    }
+
+  return retval;
+}
+
+static Octave_object
+find_nonzero_elem_idx (const Matrix& m, int nargout)
+{
+  int count = 0;
+  int m_nr = m.rows ();
+  int m_nc = m.columns ();
+
+  int i, j;
+  for (j = 0; j < m_nc; j++)
+    for (i = 0; i < m_nr; i++)
+      if (m.elem (i, j) != 0.0)
+	count++;
+
+  Matrix result;
+  Octave_object retval (nargout, result);
+
+  if (count == 0)
+    return retval;
+
+  ColumnVector i_idx (count);
+  ColumnVector j_idx (count);
+  ColumnVector v (count);
+
+  count = 0;
+  for (j = 0; j < m_nc; j++)
+    for (i = 0; i < m_nr; i++)
+      {
+	double d = m.elem (i, j);
+	if (d != 0.0)
+	  {
+	    i_idx (count) = i + 1;
+	    j_idx (count) = j + 1;
+	    v (count) = d;
+	    count++;
+	  }
+      }
+
+  tree_constant tmp (v, 1);
+  return find_to_fortran_idx (i_idx, j_idx, tmp, m_nr, m_nc, nargout);
+}
+
+static Octave_object
+find_nonzero_elem_idx (const ComplexMatrix& m, int nargout)
+{
+  int count = 0;
+  int m_nr = m.rows ();
+  int m_nc = m.columns ();
+
+  int i, j;
+  for (j = 0; j < m_nc; j++)
+    for (i = 0; i < m_nr; i++)
+      if (m.elem (i, j) != 0.0)
+	count++;
+
+  Matrix result;
+  Octave_object retval (nargout, result);
+
+  if (count == 0)
+    return retval;
+
+  ColumnVector i_idx (count);
+  ColumnVector j_idx (count);
+  ComplexColumnVector v (count);
+
+  count = 0;
+  for (j = 0; j < m_nc; j++)
+    for (i = 0; i < m_nr; i++)
+      {
+	Complex c = m.elem (i, j);
+	if (c != 0.0)
+	  {
+	    i_idx (count) = i;
+	    j_idx (count) = j;
+	    v (count) = c;
+	    count++;
+	  }
+      }
+
+  tree_constant tmp (v, 1);
+  return find_to_fortran_idx (i_idx, j_idx, tmp, m_nr, m_nc, nargout);
+}
+
+Octave_object
+find_nonzero_elem_idx (const tree_constant& a, int nargout)
+{
+  Matrix result;
+
+  nargout = (nargout == 0) ? 1 : nargout;
+  Octave_object retval (nargout, result);
+
+  tree_constant tmp = a.make_numeric ();
+
+  switch (tmp.const_type ())
+    {
+    case tree_constant_rep::matrix_constant:
+      {
+	Matrix m = tmp.matrix_value ();
+	return find_nonzero_elem_idx (m, nargout);
+      }
+      break;
+    case tree_constant_rep::scalar_constant:
+      {
+	double d = tmp.double_value ();
+	if (d != 0.0)
+	  {
+	    retval(0) = 1.0;
+	    if (nargout > 1)
+	      retval(1) = 1.0;
+	    if (nargout > 2)
+	      retval(2) = d;
+	  }
+      }
+      break;
+    case tree_constant_rep::complex_matrix_constant:
+      {
+	ComplexMatrix m = tmp.complex_matrix_value ();
+	return find_nonzero_elem_idx (m, nargout);
+      }
+      break;
+    case tree_constant_rep::complex_scalar_constant:
+      {
+	Complex c = tmp.complex_value ();
+	if (c != 0.0)
+	  {
+	    retval(0) = 1.0;
+	    if (nargout > 1)
+	      retval(1) = 1.0;
+	    if (nargout > 2)
+	      retval(2) = c;
+	  }
+      }
+      break;
+    default:
+      break;
+    }
+  return retval;
+}
+
+/*
+;;; Local Variables: ***
+;;; mode: C++ ***
+;;; page-delimiter: "^/\\*" ***
+;;; End: ***
+*/

new file mode 100644
--- /dev/null
+++ b/src/log.cc
@@ -0,0 +1,264 @@
+// f-log.cc                                           -*- C++ -*-
+/*
+
+Copyright (C) 1994 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 "EIG.h"
+
+#include "tree-const.h"
+#include "user-prefs.h"
+#include "error.h"
+#include "gripes.h"
+#include "f-log.h"
+
+// XXX FIXME XXX -- the next two functions (and expm) should really be just
+// one...
+
+Octave_object
+matrix_log (const tree_constant& a)
+{
+  Octave_object retval (1);
+
+  tree_constant tmp = a.make_numeric ();;
+    
+  if (tmp.rows () == 0 || tmp.columns () == 0)
+    {
+      int flag = user_pref.propagate_empty_matrices;
+      if (flag != 0)
+	{
+	  if (flag < 0)
+	    gripe_empty_arg ("logm", 0);
+	  Matrix m;
+	  retval(0) = m;
+	  return retval;
+	}
+      else
+	gripe_empty_arg ("logm", 1);
+    }
+
+  switch (tmp.const_type ())
+    {
+    case tree_constant_rep::matrix_constant:
+      {
+	Matrix m = tmp.matrix_value ();
+
+	int nr = m.rows ();
+	int nc = m.columns ();
+
+	if (nr == 0 || nc == 0 || nr != nc)
+	  gripe_square_matrix_required ("logm");
+	else
+	  {
+	    EIG m_eig (m);
+	    ComplexColumnVector lambda (m_eig.eigenvalues ());
+	    ComplexMatrix Q (m_eig.eigenvectors ());
+
+	    for (int i = 0; i < nr; i++)
+	      {
+		Complex elt = lambda.elem (i);
+		if (imag (elt) == 0.0 && real (elt) > 0.0)
+		  lambda.elem (i) = log (real (elt));
+		else
+		  lambda.elem (i) = log (elt);
+	      }
+
+	    ComplexDiagMatrix D (lambda);
+	    ComplexMatrix result = Q * D * Q.inverse ();
+
+	    retval(0) = result;
+	  }
+      }
+      break;
+    case tree_constant_rep::complex_matrix_constant:
+      {
+	ComplexMatrix m = tmp.complex_matrix_value ();
+
+	int nr = m.rows ();
+	int nc = m.columns ();
+
+	if (nr == 0 || nc == 0 || nr != nc)
+	  gripe_square_matrix_required ("logm");
+	else
+	  {
+	    EIG m_eig (m);
+	    ComplexColumnVector lambda (m_eig.eigenvalues ());
+	    ComplexMatrix Q (m_eig.eigenvectors ());
+
+	    for (int i = 0; i < nr; i++)
+	      {
+		Complex elt = lambda.elem (i);
+		if (imag (elt) == 0.0 && real (elt) > 0.0)
+		  lambda.elem (i) = log (real (elt));
+		else
+		  lambda.elem (i) = log (elt);
+	      }
+
+	    ComplexDiagMatrix D (lambda);
+	    ComplexMatrix result = Q * D * Q.inverse ();
+
+	    retval(0) = result;
+	  }
+      }
+      break;
+    case tree_constant_rep::scalar_constant:
+      {
+	double d = tmp.double_value ();
+	if (d > 0.0)
+	  retval(0) = log (d);
+	else
+	  {
+	    Complex dtmp (d);
+	    retval(0) = log (dtmp);
+	  }
+      }
+      break;
+    case tree_constant_rep::complex_scalar_constant:
+      {
+	Complex c = tmp.complex_value ();
+	retval(0) = log (c);
+      }
+      break;
+    default:
+      break;
+    }
+  return retval;
+}
+
+Octave_object
+matrix_sqrt (const tree_constant& a)
+{
+  Octave_object retval (1);
+
+  tree_constant tmp = a.make_numeric ();;
+    
+  if (tmp.rows () == 0 || tmp.columns () == 0)
+    {
+      int flag = user_pref.propagate_empty_matrices;
+      if (flag != 0)
+	{
+	  if (flag < 0)
+	    gripe_empty_arg ("sqrtm", 0);
+	  Matrix m;
+	  retval(0) = m;
+	  return retval;
+	}
+      else
+	gripe_empty_arg ("sqrtm", 1);
+    }
+
+  switch (tmp.const_type ())
+    {
+    case tree_constant_rep::matrix_constant:
+      {
+	Matrix m = tmp.matrix_value ();
+
+	int nr = m.rows ();
+	int nc = m.columns ();
+
+	if (nr == 0 || nc == 0 || nr != nc)
+	  gripe_square_matrix_required ("sqrtm");
+	else
+	  {
+	    EIG m_eig (m);
+	    ComplexColumnVector lambda (m_eig.eigenvalues ());
+	    ComplexMatrix Q (m_eig.eigenvectors ());
+
+	    for (int i = 0; i < nr; i++)
+	      {
+		Complex elt = lambda.elem (i);
+		if (imag (elt) == 0.0 && real (elt) > 0.0)
+		  lambda.elem (i) = sqrt (real (elt));
+		else
+		  lambda.elem (i) = sqrt (elt);
+	      }
+
+	    ComplexDiagMatrix D (lambda);
+	    ComplexMatrix result = Q * D * Q.inverse ();
+
+	    retval(0) = result;
+	  }
+      }
+      break;
+    case tree_constant_rep::complex_matrix_constant:
+      {
+	ComplexMatrix m = tmp.complex_matrix_value ();
+
+	int nr = m.rows ();
+	int nc = m.columns ();
+
+	if (nr == 0 || nc == 0 || nr != nc)
+	  gripe_square_matrix_required ("sqrtm");
+	else
+	  {
+	    EIG m_eig (m);
+	    ComplexColumnVector lambda (m_eig.eigenvalues ());
+	    ComplexMatrix Q (m_eig.eigenvectors ());
+
+	    for (int i = 0; i < nr; i++)
+	      {
+		Complex elt = lambda.elem (i);
+		if (imag (elt) == 0.0 && real (elt) > 0.0)
+		  lambda.elem (i) = sqrt (real (elt));
+		else
+		  lambda.elem (i) = sqrt (elt);
+	      }
+
+	    ComplexDiagMatrix D (lambda);
+	    ComplexMatrix result = Q * D * Q.inverse ();
+
+	    retval(0) = result;
+	  }
+      }
+      break;
+    case tree_constant_rep::scalar_constant:
+      {
+	double d = tmp.double_value ();
+	if (d > 0.0)
+	  retval(0) = sqrt (d);
+	else
+	  {
+	    Complex dtmp (d);
+	    retval(0) = sqrt (dtmp);
+	  }
+      }
+      break;
+    case tree_constant_rep::complex_scalar_constant:
+      {
+	Complex c = tmp.complex_value ();
+	retval(0) = log (c);
+      }
+      break;
+    default:
+      break;
+    }
+  return retval;
+}
+
+/*
+;;; Local Variables: ***
+;;; mode: C++ ***
+;;; page-delimiter: "^/\\*" ***
+;;; End: ***
+*/

new file mode 100644
--- /dev/null
+++ b/src/minmax.cc
@@ -0,0 +1,495 @@
+// f-minmax.cc                                           -*- C++ -*-
+/*
+
+Copyright (C) 1994 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 <math.h>
+
+#include "tree-const.h"
+#include "error.h"
+#include "f-minmax.h"
+
+#ifndef MAX
+#define MAX(a,b) ((a) > (b) ? (a) : (b))
+#endif
+
+#ifndef MIN
+#define MIN(a,b) ((a) < (b) ? (a) : (b))
+#endif
+
+static Matrix
+min (const Matrix& a, const Matrix& b)
+{
+  int nr = a.rows ();
+  int nc = a.columns ();
+  if (nr != b.rows () || nc != b.columns ())
+    {
+      error ("two-arg min expecting args of same size");
+      return Matrix ();
+    }
+
+  Matrix result (nr, nc);
+
+  for (int j = 0; j < nc; j++)
+    for (int i = 0; i < nr; i++)
+      {
+	double a_elem = a.elem (i, j);
+	double b_elem = b.elem (i, j);
+	result.elem (i, j) = MIN (a_elem, b_elem);
+      }
+
+  return result;
+}
+
+static ComplexMatrix
+min (const ComplexMatrix& a, const ComplexMatrix& b)
+{
+  int nr = a.rows ();
+  int nc = a.columns ();
+  if (nr != b.rows () || nc != b.columns ())
+    {
+      error ("two-arg min expecting args of same size");
+      return ComplexMatrix ();
+    }
+
+  ComplexMatrix result (nr, nc);
+
+  for (int j = 0; j < nc; j++)
+    for (int i = 0; i < nr; i++)
+      {
+        double abs_a_elem = abs (a.elem (i, j));
+        double abs_b_elem = abs (b.elem (i, j));
+        if (abs_a_elem < abs_b_elem)
+          result.elem (i, j) = a.elem (i, j);
+        else
+          result.elem (i, j) = b.elem (i, j);
+      }
+
+  return result;
+}
+
+static Matrix
+max (const Matrix& a, const Matrix& b)
+{
+  int nr = a.rows ();
+  int nc = a.columns ();
+  if (nr != b.rows () || nc != b.columns ())
+    {
+      error ("two-arg max expecting args of same size");
+      return Matrix ();
+    }
+
+  Matrix result (nr, nc);
+
+  for (int j = 0; j < nc; j++)
+    for (int i = 0; i < nr; i++)
+      {
+	double a_elem = a.elem (i, j);
+	double b_elem = b.elem (i, j);
+	result.elem (i, j) = MAX (a_elem, b_elem);
+      }
+
+  return result;
+}
+
+static ComplexMatrix
+max (const ComplexMatrix& a, const ComplexMatrix& b)
+{
+  int nr = a.rows ();
+  int nc = a.columns ();
+  if (nr != b.rows () || nc != b.columns ())
+    {
+      error ("two-arg max expecting args of same size");
+      return ComplexMatrix ();
+    }
+
+  ComplexMatrix result (nr, nc);
+
+  for (int j = 0; j < nc; j++)
+    for (int i = 0; i < nr; i++)
+      {
+        double abs_a_elem = abs (a.elem (i, j));
+        double abs_b_elem = abs (b.elem (i, j));
+        if (abs_a_elem > abs_b_elem)
+          result.elem (i, j) = a.elem (i, j);
+        else
+          result.elem (i, j) = b.elem (i, j);
+      }
+
+  return result;
+}
+
+Octave_object
+column_min (const Octave_object& args, int nargout)
+{
+  Octave_object retval;
+
+  tree_constant arg1;
+  tree_constant arg2;
+  tree_constant_rep::constant_type arg1_type =
+    tree_constant_rep::unknown_constant;
+  tree_constant_rep::constant_type arg2_type =
+    tree_constant_rep::unknown_constant;
+
+  int nargin = args.length ();
+
+  switch (nargin)
+    {
+    case 3:
+      arg2 = args(2).make_numeric ();
+      arg2_type = arg2.const_type ();
+// Fall through...
+    case 2:
+      arg1 = args(1).make_numeric ();
+      arg1_type = arg1.const_type ();
+      break;
+    default:
+      panic_impossible ();
+      break;
+    }
+
+  if (nargin == 2 && (nargout == 1 || nargout == 0))
+    {
+      retval.resize (1);
+      switch (arg1_type)
+	{
+        case tree_constant_rep::scalar_constant:
+	  retval(0) = arg1.double_value ();
+          break;
+        case tree_constant_rep::complex_scalar_constant:
+          retval(0) = arg1.complex_value ();
+          break;
+        case tree_constant_rep::matrix_constant:
+          {
+  	    Matrix m = arg1.matrix_value ();
+	    if (m.rows () == 1)
+	      retval(0) = m.row_min ();
+	    else
+	      retval(0) = tree_constant (m.column_min (), 0);
+ 	  }
+          break;
+        case tree_constant_rep::complex_matrix_constant:
+          {
+            ComplexMatrix m = arg1.complex_matrix_value ();
+            if (m.rows () == 1)
+              retval(0) = m.row_min ();
+            else
+              retval(0) = tree_constant (m.column_min (), 0);
+          }
+	  break;
+	default:
+	  panic_impossible ();
+	  break;
+	}
+    }
+  else if (nargin == 2 && nargout == 2)
+    {
+      retval.resize (2);
+      switch (arg1_type)
+        {
+	case tree_constant_rep::scalar_constant:
+	  {
+	    retval(0) = arg1.double_value ();
+	    retval(1) = 1;
+	  }
+          break;
+	case tree_constant_rep::complex_scalar_constant:
+	  {
+	    retval(0) = arg1.complex_value ();
+	    retval(1) = 1;
+	  }
+          break;
+	case tree_constant_rep::matrix_constant:
+	  {
+	    Matrix m = arg1.matrix_value ();
+	    if (m.rows () == 1)
+	      {
+		retval(0) = m.row_min ();
+		retval(1) = m.row_min_loc ();
+	      }
+	    else
+	      {
+		retval(0) = tree_constant (m.column_min (), 0);
+		retval(1) = tree_constant (m.column_min_loc (), 0);
+	      }
+	  }
+          break;
+	case tree_constant_rep::complex_matrix_constant:
+	  {
+	    ComplexMatrix m = arg1.complex_matrix_value ();
+	    if (m.rows () == 1)
+	      {
+		retval(0) = m.row_min ();
+		retval(1) = m.row_min_loc ();
+	      }
+	    else
+	      {
+		retval(0) = tree_constant (m.column_min (), 0);
+		retval(1) = tree_constant (m.column_min_loc (), 0);
+	      }
+	  }
+          break;
+	default:
+	  panic_impossible ();
+	  break;
+        }
+    }
+  else if (nargin == 3)
+    {
+      if (arg1.rows () == arg2.rows ()
+	  && arg1.columns () == arg2.columns ())
+	{
+	  retval.resize (1);
+          switch (arg1_type)
+            {
+	    case tree_constant_rep::scalar_constant:
+	      {
+		double result;
+		double a_elem = arg1.double_value ();
+		double b_elem = arg2.double_value ();
+		result = MIN (a_elem, b_elem);
+		retval(0) = result;
+	      }
+              break;
+	    case tree_constant_rep::complex_scalar_constant:
+	      {
+		Complex result;
+		Complex a_elem = arg1.complex_value ();
+		Complex b_elem = arg2.complex_value ();
+		if (abs (a_elem) < abs (b_elem))
+		  result = a_elem;
+		else
+		  result = b_elem;
+		retval(0) = result;
+	      }
+              break;
+	    case tree_constant_rep::matrix_constant:
+	      {
+		Matrix result;
+		result = min (arg1.matrix_value (), arg2.matrix_value ());
+		retval(0) = result;
+	      }
+              break;
+	    case tree_constant_rep::complex_matrix_constant:
+	      {
+		ComplexMatrix result;
+		result = min (arg1.complex_matrix_value (),
+			      arg2.complex_matrix_value ());
+		retval(0) = result;
+	      }
+	      break;
+	    default:
+	      panic_impossible ();
+	      break;
+	    }
+	}
+      else
+	error ("min: nonconformant matrices");
+    }
+  else
+    panic_impossible ();
+
+  return retval;
+}
+
+Octave_object
+column_max (const Octave_object& args, int nargout)
+{
+  Octave_object retval;
+
+  tree_constant arg1;
+  tree_constant arg2;
+  tree_constant_rep::constant_type arg1_type =
+    tree_constant_rep::unknown_constant;
+  tree_constant_rep::constant_type arg2_type =
+    tree_constant_rep::unknown_constant;
+
+  int nargin = args.length ();
+
+  switch (nargin)
+    {
+    case 3:
+      arg2 = args(2).make_numeric ();
+      arg2_type = arg2.const_type ();
+// Fall through...
+    case 2:
+      arg1 = args(1).make_numeric ();
+      arg1_type = arg1.const_type ();
+      break;
+    default:
+      panic_impossible ();
+      break;
+    }
+
+  if (nargin == 2 && (nargout == 1 || nargout == 0))
+    {
+      retval.resize (1);
+      switch (arg1_type)
+	{
+        case tree_constant_rep::scalar_constant:
+	  retval(0) = arg1.double_value ();
+          break;
+        case tree_constant_rep::complex_scalar_constant:
+          retval(0) = arg1.complex_value ();
+          break;
+        case tree_constant_rep::matrix_constant:
+          {
+  	    Matrix m = arg1.matrix_value ();
+	    if (m.rows () == 1)
+	      retval(0) = m.row_max ();
+	    else
+	      retval(0) = tree_constant (m.column_max (), 0);
+ 	  }
+          break;
+        case tree_constant_rep::complex_matrix_constant:
+          {
+            ComplexMatrix m = arg1.complex_matrix_value ();
+            if (m.rows () == 1)
+              retval(0) = m.row_max ();
+            else
+              retval(0) = tree_constant (m.column_max (), 0);
+          }
+	  break;
+	default:
+	  panic_impossible ();
+	  break;
+	}
+    }
+  else if (nargin == 2 && nargout == 2)
+    {
+      retval.resize (2);
+      switch (arg1_type)
+        {
+	case tree_constant_rep::scalar_constant:
+	  {
+	    retval(0) = arg1.double_value ();
+	    retval(1) = 1;
+	  }
+          break;
+	case tree_constant_rep::complex_scalar_constant:
+	  {
+	    retval(0) = arg1.complex_value ();
+	    retval(1) = 1;
+	  }
+          break;
+	case tree_constant_rep::matrix_constant:
+	  {
+	    Matrix m = arg1.matrix_value ();
+	    if (m.rows () == 1)
+	      {
+		retval(0) = m.row_max ();
+		retval(1) = m.row_max_loc ();
+	      }
+	    else
+	      {
+		retval(0) = tree_constant (m.column_max (), 0);
+		retval(1) = tree_constant (m.column_max_loc (), 0);
+	      }
+	  }
+          break;
+	case tree_constant_rep::complex_matrix_constant:
+	  {
+	    ComplexMatrix m = arg1.complex_matrix_value ();
+	    if (m.rows () == 1)
+	      {
+		retval(0) = m.row_max ();
+		retval(1) = m.row_max_loc ();
+	      }
+	    else
+	      {
+		retval(0) = tree_constant (m.column_max (), 0);
+		retval(1) = tree_constant (m.column_max_loc (), 0);
+	      }
+	  }
+          break;
+	default:
+	  panic_impossible ();
+	  break;
+        }
+    }
+  else if (nargin == 3)
+    {
+      if (arg1.rows () == arg2.rows ()
+	  && arg1.columns () == arg2.columns ())
+	{
+	  retval.resize (1);
+          switch (arg1_type)
+            {
+	    case tree_constant_rep::scalar_constant:
+	      {
+		double result;
+		double a_elem = arg1.double_value ();
+		double b_elem = arg2.double_value ();
+		result = MAX (a_elem, b_elem);
+		retval(0) = result;
+	      }
+              break;
+	    case tree_constant_rep::complex_scalar_constant:
+	      {
+		Complex result;
+		Complex a_elem = arg1.complex_value ();
+		Complex b_elem = arg2.complex_value ();
+		if (abs (a_elem) > abs (b_elem))
+		  result = a_elem;
+		else
+		  result = b_elem;
+		retval(0) = result;
+	      }
+              break;
+	    case tree_constant_rep::matrix_constant:
+	      {
+		Matrix result;
+		result = max (arg1.matrix_value (), arg2.matrix_value ());
+		retval(0) = result;
+	      }
+              break;
+	    case tree_constant_rep::complex_matrix_constant:
+	      {
+		ComplexMatrix result;
+		result = max (arg1.complex_matrix_value (),
+			      arg2.complex_matrix_value ());
+		retval(0) = result;
+	      }
+	      break;
+	    default:
+	      panic_impossible ();
+	      break;
+	    }
+	}
+      else
+	error ("max: nonconformant matrices");
+    }
+  else
+    panic_impossible ();
+
+  return retval;
+}
+
+/*
+;;; Local Variables: ***
+;;; mode: C++ ***
+;;; page-delimiter: "^/\\*" ***
+;;; End: ***
+*/

new file mode 100644
--- /dev/null
+++ b/src/sort.cc
@@ -0,0 +1,265 @@
+// f-sort.cc                                           -*- C++ -*-
+/*
+
+Copyright (C) 1994 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 "tree-const.h"
+#include "f-sort.h"
+
+static void
+mx_sort (Matrix& m, Matrix& idx, int return_idx)
+{
+  int nr = m.rows ();
+  int nc = m.columns ();
+  idx.resize (nr, nc);
+  int i, j;
+
+  if (return_idx)
+    {
+      for (j = 0; j < nc; j++)
+	for (i = 0; i < nr; i++)
+	  idx.elem (i, j) = i+1;
+    }
+
+  for (j = 0; j < nc; j++)
+    {
+      for (int gap = nr/2; gap > 0; gap /= 2)
+	for (i = gap; i < nr; i++)
+	  for (int k = i - gap;
+	       k >= 0 && m.elem (k, j) > m.elem (k+gap, j);
+	       k -= gap)
+	    {
+	      double tmp = m.elem (k, j);
+	      m.elem (k, j) = m.elem (k+gap, j);
+	      m.elem (k+gap, j) = tmp;
+
+	      if (return_idx)
+		{
+		  double tmp = idx.elem (k, j);
+		  idx.elem (k, j) = idx.elem (k+gap, j);
+		  idx.elem (k+gap, j) = tmp;
+		}
+	    }
+    }
+}
+
+static void
+mx_sort (RowVector& v, RowVector& idx, int return_idx)
+{
+  int n = v.capacity ();
+  idx.resize (n);
+  int i;
+
+  if (return_idx)
+    for (i = 0; i < n; i++)
+      idx.elem (i) = i+1;
+
+  for (int gap = n/2; gap > 0; gap /= 2)
+    for (i = gap; i < n; i++)
+      for (int k = i - gap;
+	   k >= 0 && v.elem (k) > v.elem (k+gap);
+	   k -= gap)
+	{
+	  double tmp = v.elem (k);
+	  v.elem (k) = v.elem (k+gap);
+	  v.elem (k+gap) = tmp;
+
+	  if (return_idx)
+	    {
+	      double tmp = idx.elem (k);
+	      idx.elem (k) = idx.elem (k+gap);
+	      idx.elem (k+gap) = tmp;
+	    }
+	}
+}
+
+static void
+mx_sort (ComplexMatrix& cm, Matrix& idx, int return_idx)
+{
+  int nr = cm.rows ();
+  int nc = cm.columns ();
+  idx.resize (nr, nc);
+  int i, j;
+
+  if (return_idx)
+    {
+      for (j = 0; j < nc; j++)
+	for (i = 0; i < nr; i++)
+	  idx.elem (i, j) = i+1;
+    }
+
+  for (j = 0; j < nc; j++)
+    {
+      for (int gap = nr/2; gap > 0; gap /= 2)
+	for (i = gap; i < nr; i++)
+	  for (int k = i - gap;
+	       k >= 0 && abs (cm.elem (k, j)) > abs (cm.elem (k+gap, j));
+	       k -= gap)
+	    {
+	      Complex ctmp = cm.elem (k, j);
+	      cm.elem (k, j) = cm.elem (k+gap, j);
+	      cm.elem (k+gap, j) = ctmp;
+
+	      if (return_idx)
+		{
+		  double tmp = idx.elem (k, j);
+		  idx.elem (k, j) = idx.elem (k+gap, j);
+		  idx.elem (k+gap, j) = tmp;
+		}
+	    }
+    }
+}
+
+static void
+mx_sort (ComplexRowVector& cv, RowVector& idx, int return_idx)
+{
+  int n = cv.capacity ();
+  idx.resize (n);
+  int i;
+
+  if (return_idx)
+    for (i = 0; i < n; i++)
+      idx.elem (i) = i+1;
+
+  for (int gap = n/2; gap > 0; gap /= 2)
+    for (i = gap; i < n; i++)
+      for (int k = i - gap;
+	   k >= 0 && abs (cv.elem (k)) > abs (cv.elem (k+gap));
+	   k -= gap)
+	{
+	  Complex tmp = cv.elem (k);
+	  cv.elem (k) = cv.elem (k+gap);
+	  cv.elem (k+gap) = tmp;
+
+	  if (return_idx)
+	    {
+	      double tmp = idx.elem (k);
+	      idx.elem (k) = idx.elem (k+gap);
+	      idx.elem (k+gap) = tmp;
+	    }
+	}
+}
+
+Octave_object
+sort (const Octave_object& args, int nargout)
+{
+// Assumes that we have been given the correct number of arguments.
+
+  Octave_object retval;
+
+  int return_idx = nargout > 1;
+  if (return_idx)
+    retval.resize (2);
+  else
+    retval.resize (1);
+
+  switch (args(1).const_type ())
+    {
+    case tree_constant_rep::scalar_constant:
+      {
+	retval(0) = args(1).double_value ();
+	if (return_idx)
+	  retval(1) = 1.0;
+      }
+      break;
+    case tree_constant_rep::complex_scalar_constant:
+      {
+	retval(0) = args(1).complex_value ();
+	if (return_idx)
+	  retval(1) = 1.0;
+      }
+      break;
+    case tree_constant_rep::string_constant:
+    case tree_constant_rep::range_constant:
+    case tree_constant_rep::matrix_constant:
+      {
+	Matrix m = args(1).to_matrix ();
+	if (m.rows () == 1)
+	  {
+	    int nc = m.columns ();
+	    RowVector v (nc);
+	    for (int i = 0; i < nc; i++)
+	      v.elem (i) = m.elem (0, i);
+	    RowVector idx;
+	    mx_sort (v, idx, return_idx);
+
+	    retval(0) = tree_constant (v, 0);
+	    if (return_idx)
+	      retval(1) = tree_constant (idx, 0);
+	  }
+	else
+	  {
+// Sorts m in place, optionally computes index Matrix.
+	    Matrix idx;
+	    mx_sort (m, idx, return_idx);
+
+	    retval(0) = m;
+	    if (return_idx)
+	      retval(1) = idx;
+	  }
+      }
+      break;
+    case tree_constant_rep::complex_matrix_constant:
+      {
+	ComplexMatrix cm = args(1).complex_matrix_value ();
+	if (cm.rows () == 1)
+	  {
+	    int nc = cm.columns ();
+	    ComplexRowVector cv (nc);
+	    for (int i = 0; i < nc; i++)
+	      cv.elem (i) = cm.elem (0, i);
+	    RowVector idx;
+	    mx_sort (cv, idx, return_idx);
+
+	    retval(0) = tree_constant (cv, 0);
+	    if (return_idx)
+	      retval(1) = tree_constant (idx, 0);
+	  }
+	else
+	  {
+// Sorts cm in place, optionally computes index Matrix.
+	    Matrix idx;
+	    mx_sort (cm, idx, return_idx);
+
+	    retval(0) = cm;
+	    if (return_idx)
+	      retval(1) = idx;
+	  }
+      }
+      break;
+    default:
+      panic_impossible ();
+      break;
+    }
+
+  return retval;
+}
+
+/*
+;;; Local Variables: ***
+;;; mode: C++ ***
+;;; page-delimiter: "^/\\*" ***
+;;; End: ***
+*/