view src/pt-const.cc @ 1315:611d403c7f3d

[project @ 1995-06-25 19:56:32 by jwe]
author jwe
date Sun, 25 Jun 1995 19:56:32 +0000
parents 484e94579182
children c9aecc02c6d4
line wrap: on
line source

// tree-const.cc                                         -*- C++ -*-
/*

Copyright (C) 1992, 1993, 1994, 1995 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, 59 Temple Place - Suite 330, Boston, MA  02111-1307, USA.

*/

#if defined (__GNUG__)
#pragma implementation
#endif

#ifdef HAVE_CONFIG_H
#include <config.h>
#endif

#include <string.h>
#include <iostream.h>
#include <strstream.h>

#include "tree-const.h"
#include "user-prefs.h"
#include "pager.h"
#include "error.h"
#include "gripes.h"
#include "oct-map.h"

// The following three variables could be made static members of the
// tree_constant class.

// Pointer to the blocks of memory we manage.
static tree_constant *tc_newlist = 0;

// Multiplier for allocating new blocks.
static const int tc_newlist_grow_size = 128;

Octave_map
tree_constant::map_value (void) const
{
  return rep->map_value ();
}

tree_constant::~tree_constant (void)
{
#if defined (MDEBUG)
  cerr << "~tree_constant: rep: " << rep
       << " rep->count: " << rep->count << "\n";
#endif

  if (--rep->count <= 0)
    {
      delete rep;
      rep = 0;
    }
}

void *
tree_constant::operator new (size_t size)
{
  assert (size == sizeof (tree_constant));

  if (! tc_newlist)
    {
      int block_size = tc_newlist_grow_size * sizeof (tree_constant);
      tc_newlist = (tree_constant *) new char [block_size];

      for (int i = 0; i < tc_newlist_grow_size - 1; i++)
	tc_newlist[i].freeptr = &tc_newlist[i+1];

      tc_newlist[i].freeptr = 0;
    }

  tree_constant *tmp = tc_newlist;
  tc_newlist = tc_newlist->freeptr;
  return tmp;
}

void
tree_constant::operator delete (void *p, size_t size)
{
  tree_constant *tmp = (tree_constant *) p;
  tmp->freeptr = tc_newlist;
  tc_newlist = tmp;
}

// Simple assignment.

tree_constant
tree_constant::operator = (const tree_constant& a)
{
  if (rep != a.rep)
    {
      if (--rep->count <= 0)
	delete rep;
      rep = a.rep;
      rep->count++;
    }
  return *this;  
}

tree_constant
tree_constant::lookup_map_element (const char *ref, int insert,
				   int silent)
{
  tree_constant retval;

  if (ref)
    {
      char *tmp = strsave (ref);

      SLList<char *> list;

      char *beg = tmp;
      char *end = 0;
      do
	{
	  end = strchr (beg, '.');
	  if (end)
	    *end = '\0';

	  list.append (strsave (beg));
	}
      while (end && (beg = end + 1));

      retval = lookup_map_element (list, insert, silent);

      delete [] tmp;
    }

  return retval;
}

tree_constant
tree_constant::lookup_map_element (SLList<char*>& list, int insert,
				   int silent)
{
  tree_constant retval;

  tree_constant_rep *tmp_rep = rep;

  Pix p = list.first ();
  while (p)
    {
      char *elt = list (p);

      list.next (p);

      tree_constant tmp;

      tmp = tmp_rep->lookup_map_element (elt, insert, silent);

      if (error_state)
	break;

      tmp_rep = tmp.rep;

      if (! p)
	retval = tmp;
    }

  return retval;
}

void
tree_constant::print (void)
{
  ostrstream output_buf;
  print (output_buf);
  output_buf << ends;
  maybe_page_output (output_buf);
}

// Simple structure assignment.

void
tree_constant::make_unique (void)
{
  if (rep->count > 1)
    {
      --rep->count;
      rep = new tree_constant_rep (*rep);
      rep->count = 1;
    }

  if (rep->is_map ())
    {
      for (Pix p = rep->a_map->first (); p != 0; rep->a_map->next (p))
	{
	  rep->a_map->contents (p) . make_unique ();
	}
    }
}

tree_constant::tree_constant_rep *
tree_constant::make_unique_map (void)
{
  if (! rep->is_map ())
    {
      if (--rep->count <= 0)
	delete rep;

      Octave_map m;
      rep = new tree_constant_rep (m);
      rep->count = 1;
    }

  make_unique ();

  return rep;
}

tree_constant
tree_constant::assign_map_element (SLList<char*>& list,
				   tree_constant& rhs)
{
  tree_constant_rep *tmp_rep = make_unique_map ();

  if (rhs.is_map ())
    rhs.make_unique ();

  Pix p = list.first ();
  while (p)
    {
      char *elt = list (p);

      list.next (p);

      tree_constant& tmp = tmp_rep->lookup_map_element (elt, 1);

      if (! p)
	{
	  tmp = rhs;
	  return tmp;
	}

      tmp_rep = tmp.make_unique_map ();
    }

  return tree_constant ();
}

// Indexed structure assignment.

tree_constant
tree_constant::assign_map_element (SLList<char*>& list,
				   tree_constant& rhs,
				   const Octave_object& args)
{
  tree_constant_rep *tmp_rep = make_unique_map ();

  if (rhs.is_map ())
    rhs.make_unique ();

  Pix p = list.first ();
  while (p)
    {
      char *elt = list (p);

      list.next (p);

      tree_constant& tmp = tmp_rep->lookup_map_element (elt, 1);

      if (! p)
	{
	  tmp.assign (rhs, args);
	  return tmp;
	}

      tmp_rep = tmp.make_unique_map ();
    }

  return tree_constant ();
}

void
tree_constant::print_code (ostream& os)
{
  print_code_indent (os);

  if (in_parens)
    os << "(";

  if (rep)
    rep->print_code (os);

  if (in_parens)
    os << ")";
}

int
print_as_scalar (const tree_constant& val)
{
  int nr = val.rows ();
  int nc = val.columns ();
  return (val.is_scalar_type ()
	  || val.is_string ()
	  || (val.is_matrix_type ()
	      && ((nr == 1 && nc == 1)
		  || nr == 0
		  || nc == 0)));
}

int
print_as_structure (const tree_constant& val)
{
  return val.is_map ();
}

// Construct return vector of empty matrices.  Return empty matrices
// and/or gripe when appropriate.

Octave_object
vector_of_empties (int nargout, const char *fcn_name)
{
  Octave_object retval;

// Got an empty argument, check if should gripe/return empty values.

  int flag = user_pref.propagate_empty_matrices;
  if (flag != 0)
    {
      if (flag < 0)
	gripe_empty_arg (fcn_name, 0);

      Matrix m;
      retval.resize (nargout ? nargout : 1);
      for (int i = 0; i < nargout; i++)
	retval(i) = m;
    }
  else
    gripe_empty_arg (fcn_name, 1);

  return retval;
}

// -------------------------------------------------------------------
//
// Basic stuff for the tree-constant representation class.
//
// Leave the commented #includes below to make it easy to split this
// out again, should we want to do that.
//
// -------------------------------------------------------------------

// #ifdef HAVE_CONFIG_H
// #include <config.h>
// #endif

#include <ctype.h>
// #include <string.h>
#include <fstream.h>
// #include <iostream.h>

#include "mx-base.h"
#include "Range.h"

#include "arith-ops.h"
#include "variables.h"
#include "sysdep.h"
// #include "error.h"
// #include "gripes.h"
// #include "user-prefs.h"
#include "utils.h"
#include "pr-output.h"
// #include "tree-const.h"
#include "idx-vector.h"
#include "unwind-prot.h"
// #include "oct-map.h"

#include "tc-inlines.h"

// The following three variables could be made static members of the
// TC_REP class.

// Pointer to the blocks of memory we manage.
static TC_REP *tc_rep_newlist = 0;

// Multiplier for allocating new blocks.
static const int tc_rep_newlist_grow_size = 128;

// Indentation level for structures.
static int structure_indent_level = 0;

static void
increment_structure_indent_level (void)
{
  structure_indent_level += 2;
}

static void
decrement_structure_indent_level (void)
{
  structure_indent_level -= 2;
}

static int
any_element_is_complex (const ComplexMatrix& a)
{
  int nr = a.rows ();
  int nc = a.columns ();
  for (int j = 0; j < nc; j++)
    for (int i = 0; i < nr; i++)
      if (imag (a.elem (i, j)) != 0.0)
	return 1;
  return 0;
}

// The real representation of constants.

TC_REP::tree_constant_rep (void)
{
  type_tag = unknown_constant;
  orig_text = 0;
}

TC_REP::tree_constant_rep (double d)
{
  scalar = d;
  type_tag = scalar_constant;
  orig_text = 0;
}

TC_REP::tree_constant_rep (const Matrix& m)
{
  if (m.rows () == 1 && m.columns () == 1)
    {
      scalar = m.elem (0, 0);
      type_tag = scalar_constant;
    }
  else
    {
      matrix = new Matrix (m);
      type_tag = matrix_constant;
    }
  orig_text = 0;
}

TC_REP::tree_constant_rep (const DiagMatrix& d)
{
  if (d.rows () == 1 && d.columns () == 1)
    {
      scalar = d.elem (0, 0);
      type_tag = scalar_constant;
    }
  else
    {
      matrix = new Matrix (d);
      type_tag = matrix_constant;
    }
  orig_text = 0;
}

TC_REP::tree_constant_rep (const RowVector& v, int prefer_column_vector)
{
  int len = v.capacity ();
  if (len == 1)
    {
      scalar = v.elem (0);
      type_tag = scalar_constant;
    }
  else
    {
      int pcv = (prefer_column_vector < 0)
	? user_pref.prefer_column_vectors
	  : prefer_column_vector;

      if (pcv)
	{
	  Matrix m (len, 1);
	  for (int i = 0; i < len; i++)
	    m.elem (i, 0) = v.elem (i);
	  matrix = new Matrix (m);
	  type_tag = matrix_constant;
	}
      else
	{
	  Matrix m (1, len);
	  for (int i = 0; i < len; i++)
	    m.elem (0, i) = v.elem (i);
	  matrix = new Matrix (m);
	  type_tag = matrix_constant;
	}
    }
  orig_text = 0;
}

TC_REP::tree_constant_rep (const ColumnVector& v, int prefer_column_vector)
{
  int len = v.capacity ();
  if (len == 1)
    {
      scalar = v.elem (0);
      type_tag = scalar_constant;
    }
  else
    {
      int pcv = (prefer_column_vector < 0)
	? user_pref.prefer_column_vectors
	  : prefer_column_vector;

      if (pcv)
	{
	  Matrix m (len, 1);
	  for (int i = 0; i < len; i++)
	    m.elem (i, 0) = v.elem (i);
	  matrix = new Matrix (m);
	  type_tag = matrix_constant;
	}
      else
	{
	  Matrix m (1, len);
	  for (int i = 0; i < len; i++)
	    m.elem (0, i) = v.elem (i);
	  matrix = new Matrix (m);
	  type_tag = matrix_constant;
	}
    }
  orig_text = 0;
}

TC_REP::tree_constant_rep (const Complex& c)
{
  complex_scalar = new Complex (c);
  type_tag = complex_scalar_constant;
  orig_text = 0;
}

TC_REP::tree_constant_rep (const ComplexMatrix& m)
{
  if (m.rows () == 1 && m.columns () == 1)
    {
      complex_scalar = new Complex (m.elem (0, 0));
      type_tag = complex_scalar_constant;
    }
  else
    {
      complex_matrix = new ComplexMatrix (m);
      type_tag = complex_matrix_constant;
    }
  orig_text = 0;
}

TC_REP::tree_constant_rep (const ComplexDiagMatrix& d)
{
  if (d.rows () == 1 && d.columns () == 1)
    {
      complex_scalar = new Complex (d.elem (0, 0));
      type_tag = complex_scalar_constant;
    }
  else
    {
      complex_matrix = new ComplexMatrix (d);
      type_tag = complex_matrix_constant;
    }
  orig_text = 0;
}

TC_REP::tree_constant_rep (const ComplexRowVector& v,
			   int prefer_column_vector) 
{
  int len = v.capacity ();
  if (len == 1)
    {
      complex_scalar = new Complex (v.elem (0));
      type_tag = complex_scalar_constant;
    }
  else
    {
      int pcv = (prefer_column_vector < 0)
	? user_pref.prefer_column_vectors
	  : prefer_column_vector;

      if (pcv)
	{
	  ComplexMatrix m (len, 1);
	  for (int i = 0; i < len; i++)
	    m.elem (i, 0) = v.elem (i);
	  complex_matrix = new ComplexMatrix (m);
	  type_tag = complex_matrix_constant;
	}
      else
	{
	  ComplexMatrix m (1, len);
	  for (int i = 0; i < len; i++)
	    m.elem (0, i) = v.elem (i);
	  complex_matrix = new ComplexMatrix (m);
	  type_tag = complex_matrix_constant;
	}
    }
  orig_text = 0;
}

TC_REP::tree_constant_rep (const ComplexColumnVector& v, int
			   prefer_column_vector)
{
  int len = v.capacity ();
  if (len == 1)
    {
      complex_scalar = new Complex (v.elem (0));
      type_tag = complex_scalar_constant;
    }
  else
    {
      int pcv = (prefer_column_vector < 0)
	? user_pref.prefer_column_vectors
	  : prefer_column_vector;

      if (pcv)
	{
	  ComplexMatrix m (len, 1);
	  for (int i = 0; i < len; i++)
	    m.elem (i, 0) = v.elem (i);
	  complex_matrix = new ComplexMatrix (m);
	  type_tag = complex_matrix_constant;
	}
      else
	{
	  ComplexMatrix m (1, len);
	  for (int i = 0; i < len; i++)
	    m.elem (0, i) = v.elem (i);
	  complex_matrix = new ComplexMatrix (m);
	  type_tag = complex_matrix_constant;
	}
    }
  orig_text = 0;
}

TC_REP::tree_constant_rep (const char *s)
{
  string = strsave (s);
  type_tag = string_constant;
  orig_text = 0;
}

TC_REP::tree_constant_rep (double b, double l, double i)
{
  range = new Range (b, l, i);
  int nel = range->nelem ();
  if (nel > 1)
    type_tag = range_constant;
  else
    {
      delete range;
      if (nel == 1)
	{
	  scalar = b;
	  type_tag = scalar_constant;
	}
      else if (nel == 0)
	{
	  matrix = new Matrix ();
	  type_tag = matrix_constant;
	}
      else
	{
	  type_tag = unknown_constant;
	  if (nel == -1)
	    ::error ("number of elements in range exceeds INT_MAX");
	  else
	    ::error ("invalid range");
	}
    }
  orig_text = 0;
}

TC_REP::tree_constant_rep (const Range& r)
{
  int nel = r.nelem ();
  if (nel > 1)
    {
      range = new Range (r);
      type_tag = range_constant;
    }
  else if (nel == 1)
    {
      scalar = r.base ();
      type_tag = scalar_constant;
    }
  else if (nel == 0)
    {
      matrix = new Matrix ();
      type_tag = matrix_constant;
    }
  else
    {
      type_tag = unknown_constant;
      if (nel == -1)
	::error ("number of elements in range exceeds INT_MAX");
      else
	::error ("invalid range");
    }

  orig_text = 0;
}

TC_REP::tree_constant_rep (const Octave_map& m)
{
  a_map = new Octave_map (m);
  type_tag = map_constant;
  orig_text = 0;
}

TC_REP::tree_constant_rep (TC_REP::constant_type t)
{
  assert (t == magic_colon || t == all_va_args);
  type_tag = t;
  orig_text = 0;
}

TC_REP::tree_constant_rep (const tree_constant_rep& t)
{
  type_tag = t.type_tag;

  switch (t.type_tag)
    {
    case unknown_constant:
      break;

    case scalar_constant:
      scalar = t.scalar;
      break;

    case matrix_constant:
      matrix = new Matrix (*(t.matrix));
      break;

    case string_constant:
      string = strsave (t.string);
      break;

    case complex_matrix_constant:
      complex_matrix = new ComplexMatrix (*(t.complex_matrix));
      break;

    case complex_scalar_constant:
      complex_scalar = new Complex (*(t.complex_scalar));
      break;

    case range_constant:
      range = new Range (*(t.range));
      break;

    case map_constant:
      a_map = new Octave_map (*(t.a_map));
      break;

    case magic_colon:
    case all_va_args:
      break;
    }

  orig_text = strsave (t.orig_text);
}

TC_REP::~tree_constant_rep (void)
{
  switch (type_tag)
    {
    case matrix_constant:
      delete matrix;
      break;

    case complex_scalar_constant:
      delete complex_scalar;
      break;

    case complex_matrix_constant:
      delete complex_matrix;
      break;

    case string_constant:
      delete [] string;
      break;

    case range_constant:
      delete range;
      break;

    case map_constant:
      delete a_map;
      break;

    case unknown_constant:
    case scalar_constant:
    case magic_colon:
    case all_va_args:
      break;
    }

  delete [] orig_text;
}

void *
TC_REP::operator new (size_t size)
{
  assert (size == sizeof (TC_REP));

  if (! tc_rep_newlist)
    {
      int block_size = tc_rep_newlist_grow_size * sizeof (TC_REP);
      tc_rep_newlist = (TC_REP *) new char [block_size];

      for (int i = 0; i < tc_rep_newlist_grow_size - 1; i++)
	tc_rep_newlist[i].freeptr = &tc_rep_newlist[i+1];

      tc_rep_newlist[i].freeptr = 0;
    }

  TC_REP *tmp = tc_rep_newlist;
  tc_rep_newlist = tc_rep_newlist->freeptr;
  return tmp;
}

void
TC_REP::operator delete (void *p, size_t size)
{
  TC_REP *tmp = (TC_REP *) p;
  tmp->freeptr = tc_rep_newlist;
  tc_rep_newlist = tmp;
}

int
TC_REP::rows (void) const
{
  int retval = -1;

  switch (type_tag)
    {
    case scalar_constant:
    case complex_scalar_constant:
      retval = 1;
      break;

    case string_constant:
    case range_constant:
      retval = (columns () > 0);
      break;

    case matrix_constant:
      retval = matrix->rows ();
      break;

    case complex_matrix_constant:
      retval = complex_matrix->rows ();
      break;

    default:
      break;
    }

  return retval;
}

int
TC_REP::columns (void) const
{
  int retval = -1;

  switch (type_tag)
    {
    case scalar_constant:
    case complex_scalar_constant:
      retval = 1;
      break;

    case matrix_constant:
      retval = matrix->columns ();
      break;

    case complex_matrix_constant:
      retval = complex_matrix->columns ();
      break;

    case string_constant:
      retval = strlen (string);
      break;

    case range_constant:
      retval = range->nelem ();
      break;

    default:
      break;
    }

  return retval;
}

tree_constant
TC_REP::all (void) const
{
  tree_constant retval;

  if (error_state)
    return retval;

  if (! is_numeric_type ())
    {
      tree_constant tmp = make_numeric ();

      if (error_state)
	return retval;

      return tmp.all ();
    }

  switch (type_tag)
    {
    case scalar_constant:
      {
	double status = (scalar != 0.0);
	retval = tree_constant (status);
      }
      break;

    case matrix_constant:
      {
	Matrix m = matrix->all ();
	retval = tree_constant (m);
      }
      break;

    case complex_scalar_constant:
      {
	double status = (*complex_scalar != 0.0);
	retval = tree_constant (status);
      }
      break;

    case complex_matrix_constant:
      {
	Matrix m = complex_matrix->all ();
	retval = tree_constant (m);
      }
      break;

    default:
      gripe_wrong_type_arg ("all", *this);
      break;
    }

  return retval;
}

tree_constant
TC_REP::any (void) const
{
  tree_constant retval;

  if (error_state)
    return retval;

  if (! is_numeric_type ())
    {
      tree_constant tmp = make_numeric ();

      if (error_state)
	return retval;

      return tmp.any ();
    }

  switch (type_tag)
    {
    case scalar_constant:
      {
	double status = (scalar != 0.0);
	retval = tree_constant (status);
      }
      break;

    case matrix_constant:
      {
	Matrix m = matrix->any ();
	retval = tree_constant (m);
      }
      break;

    case complex_scalar_constant:
      {
	double status = (*complex_scalar != 0.0);
	retval = tree_constant (status);
      }
      break;

    case complex_matrix_constant:
      {
	Matrix m = complex_matrix->any ();
	retval = tree_constant (m);
      }
      break;

    default:
      gripe_wrong_type_arg ("any", *this);
      break;
    }

  return retval;
}

int
TC_REP::valid_as_scalar_index (void) const
{
  return (type_tag == magic_colon
	  || (type_tag == scalar_constant 
	      && ! xisnan (scalar)
	      && NINT (scalar) == 1)
	  || (type_tag == range_constant
	      && range->nelem () == 1
	      && ! xisnan (range->base ())
	      && NINT (range->base ()) == 1));
}

int
TC_REP::valid_as_zero_index (void) const
{
  return ((type_tag == scalar_constant
	   && ! xisnan (scalar)
	   && NINT (scalar) == 0)
	  || (type_tag == matrix_constant
	      && matrix->rows () == 0
	      && matrix->columns () == 0)
	  || (type_tag == range_constant
	      && range->nelem () == 1
	      && ! xisnan (range->base ())
	      && NINT (range->base ()) == 0));
}

int
TC_REP::is_true (void) const
{
  int retval = 0;

  if (error_state)
    return retval;

  if (! is_numeric_type ())
    {
      tree_constant tmp = make_numeric ();

      if (error_state)
	return retval;

      return tmp.is_true ();
    }

  switch (type_tag)
    {
    case scalar_constant:
      retval = (scalar != 0.0);
      break;

    case matrix_constant:
      {
	Matrix m = (matrix->all ()) . all ();
	retval = (m.rows () == 1
		  && m.columns () == 1
		  && m.elem (0, 0) != 0.0);
      }
      break;

    case complex_scalar_constant:
      retval = (*complex_scalar != 0.0);
      break;

    case complex_matrix_constant:
      {
	Matrix m = (complex_matrix->all ()) . all ();
	retval = (m.rows () == 1
		  && m.columns () == 1
		  && m.elem (0, 0) != 0.0);
      }
      break;

    default:
      gripe_wrong_type_arg (0, *this);
      break;
    }

  return retval;
}

static void
warn_implicit_conversion (const char *from, const char *to)
{
  warning ("implicit conversion from %s to %s", from, to);
}

double
TC_REP::double_value (int force_string_conversion) const
{
  double retval = octave_NaN;

  switch (type_tag)
    {
    case scalar_constant:
      retval = scalar;
      break;

    case matrix_constant:
      {
	if (user_pref.do_fortran_indexing && rows () > 0 && columns () > 0)
	  retval = matrix->elem (0, 0);
	else
	  gripe_invalid_conversion ("real matrix", "real scalar");
      }
      break;

    case complex_matrix_constant:
    case complex_scalar_constant:
      {
	int flag = user_pref.ok_to_lose_imaginary_part;

	if (flag < 0)
	  warn_implicit_conversion ("complex scalar", "real scalar");

	if (flag)
	  {
	    if (type_tag == complex_scalar_constant)
	      retval = ::real (*complex_scalar);
	    else if (type_tag == complex_matrix_constant)
	      {
		if (user_pref.do_fortran_indexing
		    && rows () > 0 && columns () > 0)
		  retval = ::real (complex_matrix->elem (0, 0));
		else
		  gripe_invalid_conversion ("complex matrix", "real scalar");
	      }
	    else
	      panic_impossible ();
	  }
	else
	  gripe_invalid_conversion ("complex scalar", "real scalar");
      }
      break;

    case string_constant:
      {
	int flag = force_string_conversion;
	if (! flag)
	  flag = user_pref.implicit_str_to_num_ok;

	if (flag < 0)
	  warn_implicit_conversion ("string", "real scalar");

	int len = strlen (string);
	if (flag && (len == 1 || (len > 1 && user_pref.do_fortran_indexing)))
	  retval = toascii ((int) string[0]);
	else
	  gripe_invalid_conversion ("string", "real scalar");
      }
      break;

    case range_constant:
      {
	int nel = range->nelem ();
	if (nel == 1 || (nel > 1 && user_pref.do_fortran_indexing))
	  retval = range->base ();
	else
	  gripe_invalid_conversion ("range", "real scalar");
      }
      break;

    default:
      gripe_invalid_conversion (type_as_string (), "real scalar");
      break;
    }

  return retval;
}

Matrix
TC_REP::matrix_value (int force_string_conversion) const
{
  Matrix retval;

  switch (type_tag)
    {
    case scalar_constant:
      retval = Matrix (1, 1, scalar);
      break;

    case matrix_constant:
      retval = *matrix;
      break;

    case complex_scalar_constant:
    case complex_matrix_constant:
      {
	int flag = user_pref.ok_to_lose_imaginary_part;
	if (flag < 0)
	  warn_implicit_conversion ("complex matrix", "real matrix");

	if (flag)
	  {
	    if (type_tag == complex_scalar_constant)
	      retval = Matrix (1, 1, ::real (*complex_scalar));
	    else if (type_tag == complex_matrix_constant)
	      retval = ::real (*complex_matrix);
	    else
	      panic_impossible ();
	  }
	else
	  gripe_invalid_conversion ("complex matrix", "real matrix");
      }
      break;

    case string_constant:
      {
	int flag = force_string_conversion;
	if (! flag)
	  flag = user_pref.implicit_str_to_num_ok;

	if (flag < 0)
	  warn_implicit_conversion ("string", "real matrix");

	if (flag)
	  {
	    int len = strlen (string);

	    if (len > 0)
	      {
		retval.resize (1, len);

		for (int i = 0; i < len; i++)
		  retval.elem (0, i) = toascii ((int) string[i]);
	      }
	    else
	      retval = Matrix ();
	  }
	else
	  gripe_invalid_conversion ("string", "real matrix");
      }
      break;

    case range_constant:
      retval = range->matrix_value ();
      break;

    default:
      gripe_invalid_conversion (type_as_string (), "real matrix");
      break;
    }

  return retval;
}

Complex
TC_REP::complex_value (int force_string_conversion) const
{
  Complex retval (octave_NaN, octave_NaN);

  switch (type_tag)
    {
    case complex_scalar_constant:
      retval = *complex_scalar;
      break;

    case scalar_constant:
      retval = scalar;
      break;

    case complex_matrix_constant:
    case matrix_constant:
      {
	if (user_pref.do_fortran_indexing && rows () > 0 && columns () > 0)
	  {
	    if (type_tag == complex_matrix_constant)
	      retval = complex_matrix->elem (0, 0);
	    else
	      retval = matrix->elem (0, 0);
	  }
	else
	  gripe_invalid_conversion ("real matrix", "real scalar");
      }
      break;

    case string_constant:
      {
	int flag = force_string_conversion;
	if (! flag)
	  flag = user_pref.implicit_str_to_num_ok;

	if (flag < 0)
	  warn_implicit_conversion ("string", "complex scalar");

	int len = strlen (string);
	if (flag && (len == 1 || (len > 1 && user_pref.do_fortran_indexing)))
	  retval = toascii ((int) string[0]);
	else
	  gripe_invalid_conversion ("string", "complex scalar");
      }
      break;

    case range_constant:
      {
	int nel = range->nelem ();
	if (nel == 1 || (nel > 1 && user_pref.do_fortran_indexing))
	  retval = range->base ();
	else
	  gripe_invalid_conversion ("range", "complex scalar");
      }
      break;

    default:
      gripe_invalid_conversion (type_as_string (), "complex scalar");
      break;
    }

  return retval;
}

ComplexMatrix
TC_REP::complex_matrix_value (int force_string_conversion) const
{
  ComplexMatrix retval;

  switch (type_tag)
    {
    case scalar_constant:
      retval = ComplexMatrix (1, 1, Complex (scalar));
      break;

    case complex_scalar_constant:
      retval = ComplexMatrix (1, 1, *complex_scalar);
      break;

    case matrix_constant:
      retval = ComplexMatrix (*matrix);
      break;

    case complex_matrix_constant:
      retval = *complex_matrix;
      break;

    case string_constant:
      {
	int flag = force_string_conversion;
	if (! flag)
	  flag = user_pref.implicit_str_to_num_ok;

	if (flag < 0)
	  warn_implicit_conversion ("string", "complex matrix");

	if (flag)
	  {
	    int len = strlen (string);

	    retval.resize (1, len);

	    if (len > 1)
	      {
		for (int i = 0; i < len; i++)
		  retval.elem (0, i) = toascii ((int) string[i]);
	      }
	    else if (len == 1)
	      retval.elem (0, 0) = toascii ((int) string[0]);
	    else
	      panic_impossible ();
	  }
	else
	  gripe_invalid_conversion ("string", "real matrix");
      }
      break;

    case range_constant:
      retval = range->matrix_value ();
      break;

    default:
      gripe_invalid_conversion (type_as_string (), "complex matrix");
      break;
    }

  return retval;
}

char *
TC_REP::string_value (void) const
{
  if (type_tag == string_constant)
    return string;
  else
    {
      gripe_invalid_conversion (type_as_string (), "string");
      return 0;
    }
}

Range
TC_REP::range_value (void) const
{
  assert (type_tag == range_constant);
  return *range;
}

Octave_map
TC_REP::map_value (void) const
{
  assert (type_tag == map_constant);
  return *a_map;
}

tree_constant&
TC_REP::lookup_map_element (const char *name, int insert, int silent)
{
  static tree_constant retval;

  if (type_tag == map_constant)
    {
      Pix idx = a_map->seek (name);

      if (idx)
	return a_map->contents (idx);
      else if (insert)
	return (*a_map) [name];
      else if (! silent)
	error ("structure has no member `%s'", name);
    }
  else if (! silent)
    error ("invalid structure access attempted");

  return retval;
}

// This could be made more efficient by doing all the work here rather
// than relying on matrix_value() to do any possible type conversions.

ColumnVector
TC_REP::vector_value (int force_string_conversion,
		      int force_vector_conversion) const
{
  ColumnVector retval;

  Matrix m = matrix_value (force_string_conversion);

  if (error_state)
    return retval;

  int nr = m.rows ();
  int nc = m.columns ();
  if (nr == 1)
    {
      retval.resize (nc);
      for (int i = 0; i < nc; i++)
	retval.elem (i) = m (0, i);
    }
  else if (nc == 1)
    {
      retval.resize (nr);
      for (int i = 0; i < nr; i++)
	retval.elem (i) = m.elem (i, 0);
    }
  else if (nr > 0 && nc > 0
	   && (user_pref.do_fortran_indexing || force_vector_conversion))
    {
      retval.resize (nr * nc);
      int k = 0;
      for (int j = 0; j < nc; j++)
	for (int i = 0; i < nr; i++)
	  retval.elem (k++) = m.elem (i, j);
    }
  else
    gripe_invalid_conversion ("real matrix", "real vector");

  return retval;
}

// This could be made more efficient by doing all the work here rather
// than relying on complex_matrix_value() to do any possible type
// conversions.

ComplexColumnVector
TC_REP::complex_vector_value (int force_string_conversion,
			      int force_vector_conversion) const
{
  ComplexColumnVector retval;

  ComplexMatrix m = complex_matrix_value (force_string_conversion);

  if (error_state)
    return retval;

  int nr = m.rows ();
  int nc = m.columns ();
  if (nr == 1)
    {
      retval.resize (nc);
      for (int i = 0; i < nc; i++)
	retval.elem (i) = m (0, i);
    }
  else if (nc == 1)
    {
      retval.resize (nr);
      for (int i = 0; i < nr; i++)
	retval.elem (i) = m.elem (i, 0);
    }
  else if (nr > 0 && nc > 0
	   && (user_pref.do_fortran_indexing || force_vector_conversion))
    {
      retval.resize (nr * nc);
      int k = 0;
      for (int j = 0; j < nc; j++)
	for (int i = 0; i < nr; i++)
	  retval.elem (k++) = m.elem (i, j);
    }
  else
    gripe_invalid_conversion ("complex matrix", "complex vector");

  return retval;
}

tree_constant
TC_REP::convert_to_str (void) const
{
  tree_constant retval;

  switch (type_tag)
    {
    case complex_scalar_constant:
    case scalar_constant:
      {
	double d = double_value ();

	if (xisnan (d))
	  {
	    ::error ("invalid conversion from NaN to character");
	    return retval;
	  }
	else
	  {
	    int i = NINT (d);
// Warn about out of range conversions?
	    char s[2];
	    s[0] = (char) i;
	    s[1] = '\0';
	    retval = tree_constant (s);
	  }
      }
      break;

    case complex_matrix_constant:
    case matrix_constant:
      {
	if (rows () == 0 && columns () == 0)
	  {
	    char s = '\0';
	    retval = tree_constant (&s);
	  }
	else
	  {
	    ColumnVector v = vector_value ();
	    int len = v.length ();
	    if (len == 0)
	      {
		char s = '\0';
		retval = tree_constant (&s);
	      }
	    else
	      {
		char *s = new char [len+1];
		s[len] = '\0';
		for (int i = 0; i < len; i++)
		  {
		    double d = v.elem (i);

		    if (xisnan (d))
		      {
			::error ("invalid conversion from NaN to character");
			delete [] s;
			return retval;
		      }
		    else
		      {
			int ival = NINT (d);
// Warn about out of range conversions?
			s[i] = (char) ival;
		      }
		  }
		retval = tree_constant (s);
		delete [] s;
	      }
	  }
      }
      break;

    case range_constant:
      {
	Range r = range_value ();
	double b = r.base ();
	double incr = r.inc ();
	int nel = r.nelem ();
	char *s = new char [nel+1];
	s[nel] = '\0';
	for (int i = 0; i < nel; i++)
	  {
	    double d = b + i * incr;

	    if (xisnan (d))
	      {
		::error ("invalid conversion from NaN to character");
		delete [] s;
		return retval;
	      }
	    else
	      {
		int ival = NINT (d);
// Warn about out of range conversions?
		s[i] = (char) ival;
	      }
	  }
	retval = tree_constant (s);
	delete [] s;
      }
      break;

    case string_constant:
      retval = string;
      break;

    default:
      gripe_invalid_conversion (type_as_string (), "string");
      break;
    }

  return retval;
}

void
TC_REP::convert_to_row_or_column_vector (void)
{
  assert (type_tag == matrix_constant || type_tag == complex_matrix_constant);

  int nr = rows ();
  int nc = columns ();

  if (nr == 1 || nc == 1)
    return;

  int len = nr * nc;

  assert (len > 0);

  int new_nr = 1;
  int new_nc = 1;

  if (user_pref.prefer_column_vectors)
    new_nr = len;
  else
    new_nc = len;

  if (type_tag == matrix_constant)
    {
      Matrix *m = new Matrix (new_nr, new_nc);

      double *cop_out = matrix->fortran_vec ();

      for (int i = 0; i < len; i++)
	{
	  if (new_nr == 1)
	    m->elem (0, i) = *cop_out++;
	  else
	    m->elem (i, 0) = *cop_out++;
	}

      delete matrix;
      matrix = m;
    }
  else
    {
      ComplexMatrix *cm = new ComplexMatrix (new_nr, new_nc);

      Complex *cop_out = complex_matrix->fortran_vec ();

      for (int i = 0; i < len; i++)
	{
	  if (new_nr == 1)
	    cm->elem (0, i) = *cop_out++;
	  else
	    cm->elem (i, 0) = *cop_out++;
	}

      delete complex_matrix;
      complex_matrix = cm;
    }
}

void
TC_REP::force_numeric (int force_str_conv)
{
  switch (type_tag)
    {
    case scalar_constant:
    case matrix_constant:
    case complex_scalar_constant:
    case complex_matrix_constant:
      break;

    case string_constant:
      {
	if (! force_str_conv && ! user_pref.implicit_str_to_num_ok)
	  {
	    ::error ("failed to convert `%s' to a numeric type --", string);
	    ::error ("default conversion turned off");

	    return;
	  }

	int len = strlen (string);
	if (len > 1)
	  {
	    type_tag = matrix_constant;
	    Matrix *tm = new Matrix (1, len);
	    for (int i = 0; i < len; i++)
	      tm->elem (0, i) = toascii ((int) string[i]);
	    matrix = tm;
	  }
	else if (len == 1)
	  {
	    type_tag = scalar_constant;
	    scalar = toascii ((int) string[0]);
	  }
	else if (len == 0)
	  {
	    type_tag = matrix_constant;
	    matrix = new Matrix (0, 0);
	  }
	else
	  panic_impossible ();
      }
      break;

    case range_constant:
      {
	int len = range->nelem ();
	if (len > 1)
	  {
	    type_tag = matrix_constant;
	    Matrix *tm = new Matrix (1, len);
	    double b = range->base ();
	    double increment = range->inc ();
	    for (int i = 0; i < len; i++)
	      tm->elem (0, i) = b + i * increment;
	    matrix = tm;
	  }
	else if (len == 1)
	  {
	    type_tag = scalar_constant;
	    scalar = range->base ();
	  }
      }
      break;

    default:
      gripe_invalid_conversion (type_as_string (), "numeric type");
      break;
    }
}

tree_constant
TC_REP::make_numeric (int force_str_conv) const
{
  tree_constant retval;

  switch (type_tag)
    {
    case scalar_constant:
      retval = tree_constant (scalar);
      break;

    case matrix_constant:
      retval = tree_constant (*matrix);
      break;

    case complex_scalar_constant:
      retval = tree_constant (*complex_scalar);
      break;

    case complex_matrix_constant:
      retval = tree_constant (*complex_matrix);
      break;

    case string_constant:
      retval = tree_constant (string);
      retval.force_numeric (force_str_conv);
      break;

    case range_constant:
      retval = tree_constant (*range);
      retval.force_numeric (force_str_conv);
      break;

    default:
      gripe_invalid_conversion (type_as_string (), "numeric value");
      break;
    }

  return retval;
}

void
TC_REP::bump_value (tree_expression::type etype)
{
  switch (etype)
    {
    case tree_expression::increment:
      switch (type_tag)
	{
	case scalar_constant:
	  scalar++;
	  break;

	case matrix_constant:
	  *matrix = *matrix + 1.0;
	  break;

	case complex_scalar_constant:
	  *complex_scalar = *complex_scalar + 1.0;
	  break;

	case complex_matrix_constant:
	  *complex_matrix = *complex_matrix + 1.0;
	  break;

	case range_constant:
	  range->set_base (range->base () + 1.0);
	  range->set_limit (range->limit () + 1.0);
	  break;

	default:
	  gripe_wrong_type_arg ("operator ++", type_as_string ());
	  break;
	}
      break;

    case tree_expression::decrement:
      switch (type_tag)
	{
	case scalar_constant:
	  scalar--;
	  break;

	case matrix_constant:
	  *matrix = *matrix - 1.0;
	  break;

	case range_constant:
	  range->set_base (range->base () - 1.0);
	  range->set_limit (range->limit () - 1.0);
	  break;

	default:
	  gripe_wrong_type_arg ("operator --", type_as_string ());
	  break;
	}
      break;

    default:
      panic_impossible ();
      break;
    }
}

void
TC_REP::resize (int i, int j)
{
  switch (type_tag)
    {
    case matrix_constant:
      matrix->resize (i, j);
      break;

    case complex_matrix_constant:
      complex_matrix->resize (i, j);
      break;

    default:
      gripe_wrong_type_arg ("resize", type_as_string ());
      break;
    }
}

void
TC_REP::resize (int i, int j, double val)
{
  switch (type_tag)
    {
    case matrix_constant:
      matrix->resize (i, j, val);
      break;

    case complex_matrix_constant:
      complex_matrix->resize (i, j, val);
      break;

    default:
      gripe_wrong_type_arg ("resize", type_as_string ());
      break;
    }
}

void
TC_REP::maybe_resize (int i, int j)
{
  int nr = rows ();
  int nc = columns ();

  i++;
  j++;

  assert (i > 0 && j > 0);

  if (i > nr || j > nc)
    {
      if (user_pref.resize_on_range_error)
	resize (MAX (i, nr), MAX (j, nc), 0.0);
      else
	{
	  if (i > nr)
	    ::error ("row index = %d exceeds max row dimension = %d", i, nr);

	  if (j > nc)
	    ::error ("column index = %d exceeds max column dimension = %d",
		     j, nc);
	}
    }
}

void
TC_REP::maybe_resize (int i, force_orient f_orient)
{
  int nr = rows ();
  int nc = columns ();

  i++;

  assert (i >= 0 && (nr <= 1 || nc <= 1));

// This function never reduces the size of a vector, and all vectors
// have dimensions of at least 0x0.  If i is 0, it is either because
// a vector has been indexed with a vector of all zeros (in which case
// the index vector is empty and nothing will happen) or a vector has
// been indexed with 0 (an error which will be caught elsewhere).
  if (i == 0)
    return;

  if (nr <= 1 && nc <= 1 && i >= 1)
    {
      if (user_pref.resize_on_range_error)
	{
	  if (f_orient == row_orient)
	    resize (1, i, 0.0);
	  else if (f_orient == column_orient)
	    resize (i, 1, 0.0);
	  else if (user_pref.prefer_column_vectors)
	    resize (i, 1, 0.0);
	  else
	    resize (1, i, 0.0);
	}
      else
	::error ("matrix index = %d exceeds max dimension = %d", i, nc);
    }
  else if (nr == 1 && i > nc)
    {
      if (user_pref.resize_on_range_error)
	resize (1, i, 0.0);
      else
	::error ("matrix index = %d exceeds max dimension = %d", i, nc);
    }
  else if (nc == 1 && i > nr)
    {
      if (user_pref.resize_on_range_error)
	resize (i, 1, 0.0);
      else
	::error ("matrix index = %d exceeds max dimension = ", i, nc);
    }
}

void
TC_REP::stash_original_text (char *s)
{
  orig_text = strsave (s);
}

void
TC_REP::maybe_mutate (void)
{
  if (error_state)
    return;

  switch (type_tag)
    {
    case complex_scalar_constant:
      if (::imag (*complex_scalar) == 0.0)
	{
	  double d = ::real (*complex_scalar);
	  delete complex_scalar;
	  scalar = d;
	  type_tag = scalar_constant;
	}
      break;

    case complex_matrix_constant:
      if (! any_element_is_complex (*complex_matrix))
	{
	  Matrix *m = new Matrix (::real (*complex_matrix));
	  delete complex_matrix;
	  matrix = m;
	  type_tag = matrix_constant;
	}
      break;

    default:
      break;
    }

// Avoid calling rows() and columns() for things like magic_colon.

  int nr = 1;
  int nc = 1;
  if (type_tag == matrix_constant
      || type_tag == complex_matrix_constant
      || type_tag == range_constant)
    {
      nr = rows ();
      nc = columns ();
    }

  switch (type_tag)
    {
    case matrix_constant:
      if (nr == 1 && nc == 1)
	{
	  double d = matrix->elem (0, 0);
	  delete matrix;
	  scalar = d;
	  type_tag = scalar_constant;
	}
      break;

    case complex_matrix_constant:
      if (nr == 1 && nc == 1)
	{
	  Complex c = complex_matrix->elem (0, 0);
	  delete complex_matrix;
	  complex_scalar = new Complex (c);
	  type_tag = complex_scalar_constant;
	}
      break;

    case range_constant:
      if (nr == 1 && nc == 1)
	{
	  double d = range->base ();
	  delete range;
	  scalar = d;
	  type_tag = scalar_constant;
	}
      break;

    default:
      break;
    }
}

void
TC_REP::print (ostream& output_buf)
{
  if (error_state)
    return;

  switch (type_tag)
    {
    case scalar_constant:
      octave_print_internal (output_buf, scalar);
      break;

    case matrix_constant:
      octave_print_internal (output_buf, *matrix);
      break;

    case complex_scalar_constant:
      octave_print_internal (output_buf, *complex_scalar);
      break;

    case complex_matrix_constant:
      octave_print_internal (output_buf, *complex_matrix);
      break;

    case string_constant:
      output_buf << string << "\n";
      break;

    case range_constant:
      octave_print_internal (output_buf, *range);
      break;

    case map_constant:
      {
// XXX FIXME XXX -- would be nice to print the output in some standard
// order.  Maybe all substructures first, maybe alphabetize entries,
// etc.
	begin_unwind_frame ("TC_REP_print");

	unwind_protect_int (structure_indent_level);
	unwind_protect_int (user_pref.struct_levels_to_print);

	if (user_pref.struct_levels_to_print-- > 0)
	  {
	    output_buf << "{\n";

	    increment_structure_indent_level ();

	    for (Pix p = a_map->first (); p != 0; a_map->next (p))
	      {
		const char *key = a_map->key (p);
		tree_constant val = a_map->contents (p);

		output_buf.form ("%*s%s = ", structure_indent_level,
				 "", key);

		if (! (print_as_scalar (val) || print_as_structure (val))) 
		  output_buf << "\n";

		val.print (output_buf);
	      }

	    decrement_structure_indent_level ();

	    output_buf.form ("%*s%s", structure_indent_level, "", "}\n");
	  }
	else
	  output_buf << "<structure>\n";

	run_unwind_frame ("TC_REP_print");
      }
      break;

    case unknown_constant:
    case magic_colon:
    case all_va_args:
      panic_impossible ();
      break;
    }
}

void
TC_REP::print_code (ostream& os)
{
  switch (type_tag)
    {
    case scalar_constant:
      if (orig_text)
	os << orig_text;
      else
	octave_print_internal (os, scalar, 1);
      break;

    case matrix_constant:
      octave_print_internal (os, *matrix, 1);
      break;

    case complex_scalar_constant:
     {
	double re = complex_scalar->real ();
	double im = complex_scalar->imag ();

// If we have the original text and a pure imaginary, just print the
// original text, because this must be a constant that was parsed as
// part of a function.

	if (orig_text && re == 0.0 && im > 0.0)
	  os << orig_text;
	else
	  octave_print_internal (os, *complex_scalar, 1);
      }
      break;

    case complex_matrix_constant:
      octave_print_internal (os, *complex_matrix, 1);
      break;

    case string_constant:
      {
	os << "\"";
	char *s, *t = string;
	while ((s = undo_string_escape (*t++)))
	  os << s;
	os << "\"";
      }
      break;

    case range_constant:
      octave_print_internal (os, *range, 1);
      break;

    case magic_colon:
      os << ":";
      break;

    case all_va_args:
      os << "all_va_args";
      break;

    case map_constant:
    case unknown_constant:
      panic_impossible ();
      break;
    }
}

void
TC_REP::gripe_wrong_type_arg (const char *name,
			      const tree_constant_rep& tcr) const
{
  if (name)
    ::error ("%s: wrong type argument `%s'", name, tcr.type_as_string ());
  else
    ::error ("wrong type argument `%s'", name, tcr.type_as_string ());
}

char *
TC_REP::type_as_string (void) const
{
  switch (type_tag)
    {
    case scalar_constant:
      return "real scalar";

    case matrix_constant:
      return "real matrix";

    case complex_scalar_constant:
      return "complex scalar";

    case complex_matrix_constant:
      return "complex matrix";

    case string_constant:
      return "string";

    case range_constant:
      return "range";

    case map_constant:
      return "structure";

    default:
      return "<unknown type>";
    }
}

tree_constant
do_binary_op (tree_constant& a, tree_constant& b, tree_expression::type t)
{
  tree_constant retval;

  int first_empty = (a.rows () == 0 || a.columns () == 0);
  int second_empty = (b.rows () == 0 || b.columns () == 0);

  if (first_empty || second_empty)
    {
      int flag = user_pref.propagate_empty_matrices;
      if (flag < 0)
	warning ("binary operation on empty matrix");
      else if (flag == 0)
	{
	  ::error ("invalid binary operation on empty matrix");
	  return retval;
	}
    }

  tree_constant tmp_a = a.make_numeric ();

  if (error_state)
    return retval;

  tree_constant tmp_b = b.make_numeric ();

  if (error_state)
    return retval;

  TC_REP::constant_type a_type = tmp_a.const_type ();
  TC_REP::constant_type b_type = tmp_b.const_type ();

  double d1, d2;
  Matrix m1, m2;
  Complex c1, c2;
  ComplexMatrix cm1, cm2;

  switch (a_type)
    {
    case TC_REP::scalar_constant:

      d1 = tmp_a.double_value ();

      switch (b_type)
	{
	case TC_REP::scalar_constant:
	  d2 = tmp_b.double_value ();
	  retval = do_binary_op (d1, d2, t);
	  break;

	case TC_REP::matrix_constant:
	  m2 = tmp_b.matrix_value ();
	  retval = do_binary_op (d1, m2, t);
	  break;

	case TC_REP::complex_scalar_constant:
	  c2 = tmp_b.complex_value ();
	  retval = do_binary_op (d1, c2, t);
	  break;

	case TC_REP::complex_matrix_constant:
	  cm2 = tmp_b.complex_matrix_value ();
	  retval = do_binary_op (d1, cm2, t);
	  break;

	default:
	  gripe_wrong_type_arg_for_binary_op (tmp_b);
	  break;
	}
      break;

    case TC_REP::matrix_constant:

      m1 = tmp_a.matrix_value ();

      switch (b_type)
	{
	case TC_REP::scalar_constant:
	  d2 = tmp_b.double_value ();
	  retval = do_binary_op (m1, d2, t);
	  break;

	case TC_REP::matrix_constant:
	  m2 = tmp_b.matrix_value ();
	  retval = do_binary_op (m1, m2, t);
	  break;

	case TC_REP::complex_scalar_constant:
	  c2 = tmp_b.complex_value ();
	  retval = do_binary_op (m1, c2, t);
	  break;

	case TC_REP::complex_matrix_constant:
	  cm2 = tmp_b.complex_matrix_value ();
	  retval = do_binary_op (m1, cm2, t);
	  break;

	default:
	  gripe_wrong_type_arg_for_binary_op (tmp_b);
	  break;
	}
      break;

    case TC_REP::complex_scalar_constant:

      c1 = tmp_a.complex_value ();

      switch (b_type)
	{
	case TC_REP::scalar_constant:
	  d2 = tmp_b.double_value ();
	  retval = do_binary_op (c1, d2, t);
	  break;

	case TC_REP::matrix_constant:
	  m2 = tmp_b.matrix_value ();
	  retval = do_binary_op (c1, m2, t);
	  break;

	case TC_REP::complex_scalar_constant:
	  c2 = tmp_b.complex_value ();
	  retval = do_binary_op (c1, c2, t);
	  break;

	case TC_REP::complex_matrix_constant:
	  cm2 = tmp_b.complex_matrix_value ();
	  retval = do_binary_op (c1, cm2, t);
	  break;

	default:
	  gripe_wrong_type_arg_for_binary_op (tmp_b);
	  break;
	}
      break;

    case TC_REP::complex_matrix_constant:

      cm1 = tmp_a.complex_matrix_value ();

      switch (b_type)
	{
	case TC_REP::scalar_constant:
	  d2 = tmp_b.double_value ();
	  retval = do_binary_op (cm1, d2, t);
	  break;

	case TC_REP::matrix_constant:
	  m2 = tmp_b.matrix_value ();
	  retval = do_binary_op (cm1, m2, t);
	  break;

	case TC_REP::complex_scalar_constant:
	  c2 = tmp_b.complex_value ();
	  retval = do_binary_op (cm1, c2, t);
	  break;

	case TC_REP::complex_matrix_constant:
	  cm2 = tmp_b.complex_matrix_value ();
	  retval = do_binary_op (cm1, cm2, t);
	  break;

	default:
	  gripe_wrong_type_arg_for_binary_op (tmp_b);
	  break;
	}
      break;

    default:
      gripe_wrong_type_arg_for_binary_op (tmp_a);
      break;
    }

  return retval;
}

tree_constant
do_unary_op (tree_constant& a, tree_expression::type t)
{
  tree_constant retval;

  if (a.rows () == 0 || a.columns () == 0)
    {
      int flag = user_pref.propagate_empty_matrices;
      if (flag < 0)
	warning ("unary operation on empty matrix");
      else if (flag == 0)
	{
	  ::error ("invalid unary operation on empty matrix");
	  return retval;
	}
    }

  tree_constant tmp_a = a.make_numeric ();

  if (error_state)
    return retval;

  switch (tmp_a.const_type ())
    {
    case TC_REP::scalar_constant:
      retval = do_unary_op (tmp_a.double_value (), t);
      break;

    case TC_REP::matrix_constant:
      {
	Matrix m = tmp_a.matrix_value ();
	retval = do_unary_op (m, t);
      }
      break;

    case TC_REP::complex_scalar_constant:
      retval = do_unary_op (tmp_a.complex_value (), t);
      break;

    case TC_REP::complex_matrix_constant:
      {
	ComplexMatrix m = tmp_a.complex_matrix_value ();
	retval = do_unary_op (m, t);
      }
      break;

    default:
      gripe_wrong_type_arg_for_unary_op (tmp_a);
      break;
    }

  return retval;
}

// -------------------------------------------------------------------
//
// Indexing operations for the tree-constant representation class.
//
// Leave the commented #includes below to make it easy to split this
// out again, should we want to do that.
//
// -------------------------------------------------------------------

// #ifdef HAVE_CONFIG_H
// #include <config.h>
// #endif

// #include <ctype.h>
// #include <string.h>
// #include <fstream.h>
// #include <iostream.h>
// #include <strstream.h>

// #include "mx-base.h"
// #include "Range.h"

// #include "arith-ops.h"
// #include "variables.h"
// #include "sysdep.h"
// #include "error.h"
// #include "gripes.h"
// #include "user-prefs.h"
// #include "utils.h"
// #include "pager.h"
// #include "pr-output.h"
// #include "tree-const.h"
// #include "idx-vector.h"
// #include "oct-map.h"

// #include "tc-inlines.h"

// Indexing functions.

// This is the top-level indexing function.

tree_constant
TC_REP::do_index (const Octave_object& args)
{
  tree_constant retval;

  if (error_state)
    return retval;

  if (rows () == 0 || columns () == 0)
    {
      switch (args.length ())
	{
	case 2:
	  if (! args(1).is_magic_colon ()
	      && args(1).rows () != 0 && args(1).columns () != 0)
	    goto index_error;

	case 1:
	  if (! args(0).is_magic_colon ()
	      && args(0).rows () != 0 && args(0).columns () != 0)
	    goto index_error;

	  return Matrix ();

	default:
	index_error:
	  ::error ("attempt to index empty matrix");
	  return retval;
	}
    }

  switch (type_tag)
    {
    case complex_scalar_constant:
    case scalar_constant:
      retval = do_scalar_index (args);
      break;

    case complex_matrix_constant:
    case matrix_constant:
      retval = do_matrix_index (args);
      break;

    case string_constant:
      gripe_string_invalid ();
//      retval = do_string_index (args);
      break;

    default:

// This isn\'t great, but it\'s easier than implementing a lot of
// other special indexing functions.

      force_numeric ();

      if (! error_state && is_numeric_type ())
	retval = do_index (args);

      break;
    }

  return retval;
}

tree_constant
TC_REP::do_scalar_index (const Octave_object& args) const
{
  tree_constant retval;

  if (valid_scalar_indices (args))
    {
      if (type_tag == scalar_constant)
	retval = scalar;
      else if (type_tag == complex_scalar_constant)
	retval = *complex_scalar;
      else
	panic_impossible ();

      return retval;
    }
  else
    {
      int rows = -1;
      int cols = -1;

      int nargin = args.length ();

      switch (nargin)
	{
	case 2:
	  {
	    tree_constant arg = args(1);

	    if (arg.is_matrix_type ())
	      {
		Matrix mj = arg.matrix_value ();

		idx_vector j (mj, user_pref.do_fortran_indexing, "", 1);
		if (! j)
		  return retval;

		int jmax = j.max ();
		int len = j.length ();
		if (len == j.ones_count ())
		  cols = len;
		else if (jmax > 0)
		  {
		    error ("invalid scalar index = %d", jmax+1);
		    return retval;
		  }
	      }
	    else if (arg.const_type () == magic_colon)
	      {
		cols = 1;
	      }
	    else if (arg.is_scalar_type ())
	      {
		double dval = arg.double_value ();
		if (! xisnan (dval))
		  {
		    int ival = NINT (dval);
		    if (ival == 1)
		      cols = 1;
		    else if (ival == 0)
		      cols = 0;
		    else
		      break;;
		  }
		else
		  break;
	      }
	    else
	      break;
	  }

// Fall through...

	case 1:
	  {
	    tree_constant arg = args(0);

	    if (arg.is_matrix_type ())
	      {
		Matrix mi = arg.matrix_value ();

		idx_vector i (mi, user_pref.do_fortran_indexing, "", 1);
		if (! i)
		  return retval;

		int imax = i.max ();
		int len = i.length ();
		if (len == i.ones_count ())
		  rows = len;
		else if (imax > 0)
		  {
		    error ("invalid scalar index = %d", imax+1);
		    return retval;
		  }
	      }
	    else if (arg.const_type () == magic_colon)
	      {
		rows = 1;
	      }
	    else if (arg.is_scalar_type ())
	      {
		double dval = arg.double_value ();

		if (! xisnan (dval))
		  {
		    int ival = NINT (dval);
		    if (ival == 1)
		      rows = 1;
		    else if (ival == 0)
		      rows = 0;
		    else
		      break;
		  }
		else
		  break;
	      }
	    else
	      break;

// If only one index, cols will not be set, so we set it.
// If single index is [], rows will be zero, and we should set cols to
// zero too.

	    if (cols < 0)
	      {
		if (rows == 0)
		  cols = 0;
		else
		  {
		    if (user_pref.prefer_column_vectors)
		      cols = 1;
		    else
		      {
			cols = rows;
			rows = 1;
		      }
		  }
	      }

	    if (type_tag == scalar_constant)
	      {
		return Matrix (rows, cols, scalar);
	      }
	    else if (type_tag == complex_scalar_constant)
	      {
		return ComplexMatrix (rows, cols, *complex_scalar);
	      }
	    else
	      panic_impossible ();
	  }
	  break;

	default:
	  ::error ("invalid number of arguments for scalar type");
	  return tree_constant ();
	  break;
	}
    }

  ::error ("index invalid or out of range for scalar type");
  return tree_constant ();
}

tree_constant
TC_REP::do_matrix_index (const Octave_object& args) const
{
  tree_constant retval;

  int nargin = args.length ();

  switch (nargin)
    {
    case 1:
      {
	tree_constant arg = args(0);

	if (arg.is_undefined ())
	  ::error ("matrix index is a null expression");
	else
	  retval = do_matrix_index (arg);
      }
      break;

    case 2:
      {
	tree_constant arg_a = args(0);
	tree_constant arg_b = args(1);

	if (arg_a.is_undefined ())
	::error ("first matrix index is a null expression");
	else if (arg_b.is_undefined ())
	  ::error ("second matrix index is a null expression");
	else
	  retval = do_matrix_index (arg_a, arg_b);
      }
      break;

    default:
      if (nargin == 0)
	::error ("matrix indices expected, but none provided");
      else
	::error ("too many indices for matrix expression");
      break;
    }

  return  retval;
}

tree_constant
TC_REP::do_matrix_index (const tree_constant& i_arg) const
{
  tree_constant retval;

  int nr = rows ();
  int nc = columns ();

  if (user_pref.do_fortran_indexing)
    retval = fortran_style_matrix_index (i_arg);
  else if (nr <= 1 || nc <= 1)
    retval = do_vector_index (i_arg);
  else
    ::error ("single index only valid for row or column vector");

  return retval;
}

tree_constant
TC_REP::do_matrix_index (const tree_constant& i_arg,
			 const tree_constant& j_arg) const
{
  tree_constant retval;

  tree_constant tmp_i = i_arg.make_numeric_or_range_or_magic ();

  if (error_state)
    return retval;

  TC_REP::constant_type itype = tmp_i.const_type ();

  switch (itype)
    {
    case complex_scalar_constant:
    case scalar_constant:
      {
        int i = tree_to_mat_idx (tmp_i.double_value ());
	retval = do_matrix_index (i, j_arg);
      }
      break;

    case complex_matrix_constant:
    case matrix_constant:
      {
	Matrix mi = tmp_i.matrix_value ();
	idx_vector iv (mi, user_pref.do_fortran_indexing, "row", rows ());
	if (! iv)
	  return tree_constant ();

	if (iv.length () == 0)
	  {
	    Matrix mtmp;
	    retval = tree_constant (mtmp);
	  }
	else
	  retval = do_matrix_index (iv, j_arg);
      }
      break;

    case string_constant:
      gripe_string_invalid ();
      break;

    case range_constant:
      {
	Range ri = tmp_i.range_value ();
	int nr = rows ();
	if (nr == 2 && is_zero_one (ri))
	  {
	    retval = do_matrix_index (1, j_arg);
	  }
	else if (nr == 2 && is_one_zero (ri))
	  {
	    retval = do_matrix_index (0, j_arg);
	  }
	else
	  {
	    if (index_check (ri, "row") < 0)
	      return tree_constant ();
	    retval = do_matrix_index (ri, j_arg);
	  }
      }
      break;

    case magic_colon:
      retval = do_matrix_index (magic_colon, j_arg);
      break;

    default:
      panic_impossible ();
      break;
    }

  return retval;
}

tree_constant
TC_REP::do_matrix_index (TC_REP::constant_type mci) const
{
  assert (mci == magic_colon);

  tree_constant retval;
  int nr =  rows ();
  int nc =  columns ();
  int size = nr * nc;
  if (size > 0)
    {
      CRMATRIX (m, cm, size, 1);
      int idx = 0;
      for (int j = 0; j < nc; j++)
	for (int i = 0; i < nr; i++)
	  {
	    CRMATRIX_ASSIGN_REP_ELEM (m, cm, idx, 0, i, j);
	    idx++;
	  }
      ASSIGN_CRMATRIX_TO (retval, m, cm);
    }
  return retval;
}

tree_constant
TC_REP::fortran_style_matrix_index (const tree_constant& i_arg) const
{
  tree_constant retval;

  tree_constant tmp_i = i_arg.make_numeric_or_magic ();

  if (error_state)
    return retval;

  TC_REP::constant_type itype = tmp_i.const_type ();

  int nr = rows ();
  int nc = columns ();

  switch (itype)
    {
    case complex_scalar_constant:
    case scalar_constant:
      {
	double dval = tmp_i.double_value ();

	if (xisnan (dval))
	  {
	    ::error ("NaN is invalid as a matrix index");
	    return tree_constant ();
	  }
	else
	  {
	    int i = NINT (dval);
	    int ii = fortran_row (i, nr) - 1;
	    int jj = fortran_column (i, nr) - 1;
	    if (index_check (i-1, "") < 0)
	      return tree_constant ();
	    if (range_max_check (i-1, nr * nc) < 0)
	      return tree_constant ();
	    retval = do_matrix_index (ii, jj);
	  }
      }
      break;

    case complex_matrix_constant:
    case matrix_constant:
      {
	Matrix mi = tmp_i.matrix_value ();
	if (mi.rows () == 0 || mi.columns () == 0)
	  {
	    Matrix mtmp;
	    retval = tree_constant (mtmp);
	  }
	else
	  {
// Yes, we really do want to call this with mi.
	    retval = fortran_style_matrix_index (mi);
	  }
      }
      break;

    case string_constant:
      gripe_string_invalid ();
      break;

    case range_constant:
      gripe_range_invalid ();
      break;

    case magic_colon:
      retval = do_matrix_index (magic_colon);
      break;

    default:
      panic_impossible ();
      break;
    }

  return retval;
}

tree_constant
TC_REP::fortran_style_matrix_index (const Matrix& mi) const
{
  assert (is_matrix_type ());

  tree_constant retval;

  int nr = rows ();
  int nc = columns ();

  int len = nr * nc;

  int index_nr = mi.rows ();
  int index_nc = mi.columns ();

  if (index_nr >= 1 && index_nc >= 1)
    {
      const double *cop_out = 0;
      const Complex *c_cop_out = 0;
      int real_type = type_tag == matrix_constant;
      if (real_type)
	cop_out = matrix->data ();
      else
	c_cop_out = complex_matrix->data ();

      const double *cop_out_index = mi.data ();

      idx_vector iv (mi, 1, "", len);
      if (! iv || range_max_check (iv.max (), len) < 0)
	return retval;

      int result_size = iv.length ();

// XXX FIXME XXX -- there is way too much duplicate code here...

      if (iv.one_zero_only ())
	{
	  if (iv.ones_count () == 0)
	    {
	      retval = Matrix ();
	    }
	  else
	    {
	      if (nr == 1)
		{
		  CRMATRIX (m, cm, 1, result_size);

		  for (int i = 0; i < result_size; i++)
		    {
		      int idx = iv.elem (i);
		      CRMATRIX_ASSIGN_ELEM (m, cm, 0, i, cop_out [idx],
					    c_cop_out [idx], real_type);
		    }

		  ASSIGN_CRMATRIX_TO (retval, m, cm);
		}
	      else
		{
		  CRMATRIX (m, cm, result_size, 1);

		  for (int i = 0; i < result_size; i++)
		    {
		      int idx = iv.elem (i);
		      CRMATRIX_ASSIGN_ELEM (m, cm, i, 0, cop_out [idx],
					    c_cop_out [idx], real_type);
		    }

		  ASSIGN_CRMATRIX_TO (retval, m, cm);
		}
	    }
	}
      else if (nc == 1)
	{
	  CRMATRIX (m, cm, result_size, 1);

	  for (int i = 0; i < result_size; i++)
	    {
	      int idx = iv.elem (i);
	      CRMATRIX_ASSIGN_ELEM (m, cm, i, 0, cop_out [idx],
				    c_cop_out [idx], real_type);
	    }

	  ASSIGN_CRMATRIX_TO (retval, m, cm);
	}
      else if (nr == 1)
	{
	  CRMATRIX (m, cm, 1, result_size);

	  for (int i = 0; i < result_size; i++)
	    {
	      int idx = iv.elem (i);
	      CRMATRIX_ASSIGN_ELEM (m, cm, 0, i, cop_out [idx],
				    c_cop_out [idx], real_type);
	    }

	  ASSIGN_CRMATRIX_TO (retval, m, cm);
	}
      else
	{
	  CRMATRIX (m, cm, index_nr, index_nc);

	  for (int j = 0; j < index_nc; j++)
	    for (int i = 0; i < index_nr; i++)
	      {
		double tmp = *cop_out_index++;
		int idx = tree_to_mat_idx (tmp);
		CRMATRIX_ASSIGN_ELEM (m, cm, i, j, cop_out [idx],
				      c_cop_out [idx], real_type);
	      }

	  ASSIGN_CRMATRIX_TO (retval, m, cm);
	}
    }
  else
    {
      if (index_nr == 0 || index_nc == 0)
	::error ("empty matrix invalid as index");
      else
	::error ("invalid matrix index");
      return tree_constant ();
    }

  return retval;
}

tree_constant
TC_REP::do_vector_index (const tree_constant& i_arg) const
{
  tree_constant retval;

  tree_constant tmp_i = i_arg.make_numeric_or_range_or_magic ();

  if (error_state)
    return retval;

  TC_REP::constant_type itype = tmp_i.const_type ();

  int nr = rows ();
  int nc = columns ();

  int len = MAX (nr, nc);

  assert ((nr == 1 || nc == 1) && ! user_pref.do_fortran_indexing);

  int swap_indices = (nr == 1);

  switch (itype)
    {
    case complex_scalar_constant:
    case scalar_constant:
      {
        int i = tree_to_mat_idx (tmp_i.double_value ());
        if (index_check (i, "") < 0)
	  return tree_constant ();
        if (swap_indices)
          {
	    if (range_max_check (i, nc) < 0)
	      return tree_constant ();
	    retval = do_matrix_index (0, i);
          }
        else
          {
	    if (range_max_check (i, nr) < 0)
	      return tree_constant ();
	    retval = do_matrix_index (i, 0);
          }
      }
      break;

    case complex_matrix_constant:
    case matrix_constant:
      {
        Matrix mi = tmp_i.matrix_value ();
	if (mi.rows () == 0 || mi.columns () == 0)
	  {
	    Matrix mtmp;
	    retval = tree_constant (mtmp);
	  }
	else
	  {
	    idx_vector iv (mi, user_pref.do_fortran_indexing, "", len);
	    if (! iv)
	      return tree_constant ();

	    if (swap_indices)
	      {
		if (range_max_check (iv.max (), nc) < 0)
		  return tree_constant ();
		retval = do_matrix_index (0, iv);
	      }
	    else
	      {
		if (range_max_check (iv.max (), nr) < 0)
		  return tree_constant ();
		retval = do_matrix_index (iv, 0);
	      }
	  }
      }
      break;

    case string_constant:
      gripe_string_invalid ();
      break;

    case range_constant:
      {
        Range ri = tmp_i.range_value ();
	if (len == 2 && is_zero_one (ri))
	  {
	    if (swap_indices)
	      retval = do_matrix_index (0, 1);
	    else
	      retval = do_matrix_index (1, 0);
	  }
	else if (len == 2 && is_one_zero (ri))
	  {
	    retval = do_matrix_index (0, 0);
	  }
	else
	  {
	    if (index_check (ri, "") < 0)
	      return tree_constant ();
	    if (swap_indices)
	      {
		if (range_max_check (tree_to_mat_idx (ri.max ()), nc) < 0)
		  return tree_constant ();
		retval = do_matrix_index (0, ri);
	      }
	    else
	      {
		if (range_max_check (tree_to_mat_idx (ri.max ()), nr) < 0)
		  return tree_constant ();
		retval = do_matrix_index (ri, 0);
	      }
	  }
      }
      break;

    case magic_colon:
      if (swap_indices)
        retval = do_matrix_index (0, magic_colon);
      else
        retval = do_matrix_index (magic_colon, 0);
      break;

    default:
      panic_impossible ();
      break;
    }

  return retval;
}

tree_constant
TC_REP::do_matrix_index (int i, const tree_constant& j_arg) const
{
  tree_constant retval;

  tree_constant tmp_j = j_arg.make_numeric_or_range_or_magic ();

  if (error_state)
    return retval;

  TC_REP::constant_type jtype = tmp_j.const_type ();

  int nr = rows ();
  int nc = columns ();

  switch (jtype)
    {
    case complex_scalar_constant:
    case scalar_constant:
      {
	if (index_check (i, "row") < 0)
	  return tree_constant ();
	int j = tree_to_mat_idx (tmp_j.double_value ());
	if (index_check (j, "column") < 0)
	  return tree_constant ();
	if (range_max_check (i, j, nr, nc) < 0)
	  return tree_constant ();
	retval = do_matrix_index (i, j);
      }
      break;

    case complex_matrix_constant:
    case matrix_constant:
      {
	if (index_check (i, "row") < 0)
	  return tree_constant ();
	Matrix mj = tmp_j.matrix_value ();
	idx_vector jv (mj, user_pref.do_fortran_indexing, "column", nc);
	if (! jv)
	  return tree_constant ();

	if (jv.length () == 0)
	  {
	    Matrix mtmp;
	    retval = tree_constant (mtmp);
	  }
	else
	  {
	    if (range_max_check (i, jv.max (), nr, nc) < 0)
	      return tree_constant ();
	    retval = do_matrix_index (i, jv);
	  }
      }
      break;

    case string_constant:
      gripe_string_invalid ();
      break;

    case range_constant:
      {
	if (index_check (i, "row") < 0)
	  return tree_constant ();
	Range rj = tmp_j.range_value ();
	if (nc == 2 && is_zero_one (rj))
	  {
	    retval = do_matrix_index (i, 1);
	  }
	else if (nc == 2 && is_one_zero (rj))
	  {
	    retval = do_matrix_index (i, 0);
	  }
	else
	  {
	    if (index_check (rj, "column") < 0)
	      return tree_constant ();
	    if (range_max_check (i, tree_to_mat_idx (rj.max ()), nr, nc) < 0)
	      return tree_constant ();
	    retval = do_matrix_index (i, rj);
	  }
      }
      break;

    case magic_colon:
      if (i == -1 && nr == 1)
	return Matrix ();
      if (index_check (i, "row") < 0
	  || range_max_check (i, 0, nr, nc) < 0)
	return tree_constant ();
      retval = do_matrix_index (i, magic_colon);
      break;

    default:
      panic_impossible ();
      break;
    }

  return retval;
}

tree_constant
TC_REP::do_matrix_index (const idx_vector& iv,
			 const tree_constant& j_arg) const
{
  tree_constant retval;

  tree_constant tmp_j = j_arg.make_numeric_or_range_or_magic ();

  if (error_state)
    return retval;

  TC_REP::constant_type jtype = tmp_j.const_type ();

  int nr = rows ();
  int nc = columns ();

  switch (jtype)
    {
    case complex_scalar_constant:
    case scalar_constant:
      {
	int j = tree_to_mat_idx (tmp_j.double_value ());
	if (index_check (j, "column") < 0)
	  return tree_constant ();
	if (range_max_check (iv.max (), j, nr, nc) < 0)
	  return tree_constant ();
	retval = do_matrix_index (iv, j);
      }
      break;

    case complex_matrix_constant:
    case matrix_constant:
      {
	Matrix mj = tmp_j.matrix_value ();
	idx_vector jv (mj, user_pref.do_fortran_indexing, "column", nc);
	if (! jv)
	  return tree_constant ();

	if (jv.length () == 0)
	  {
	    Matrix mtmp;
	    retval = tree_constant (mtmp);
	  }
	else
	  {
	    if (range_max_check (iv.max (), jv.max (), nr, nc) < 0)
	      return tree_constant ();
	    retval = do_matrix_index (iv, jv);
	  }
      }
      break;

    case string_constant:
      gripe_string_invalid ();
      break;

    case range_constant:
      {
	Range rj = tmp_j.range_value ();
	if (nc == 2 && is_zero_one (rj))
	  {
	    retval = do_matrix_index (iv, 1);
	  }
	else if (nc == 2 && is_one_zero (rj))
	  {
	    retval = do_matrix_index (iv, 0);
	  }
	else
	  {
	    if (index_check (rj, "column") < 0)
	      return tree_constant ();
	    if (range_max_check (iv.max (), tree_to_mat_idx (rj.max ()),
				 nr, nc) < 0)
	      return tree_constant ();
	    retval = do_matrix_index (iv, rj);
	  }
      }
      break;

    case magic_colon:
      if (range_max_check (iv.max (), 0, nr, nc) < 0)
	return tree_constant ();
      retval = do_matrix_index (iv, magic_colon);
      break;

    default:
      panic_impossible ();
      break;
    }

  return retval;
}

tree_constant
TC_REP::do_matrix_index (const Range& ri,
			 const tree_constant& j_arg) const
{
  tree_constant retval;

  tree_constant tmp_j = j_arg.make_numeric_or_range_or_magic ();

  if (error_state)
    return retval;

  TC_REP::constant_type jtype = tmp_j.const_type ();

  int nr = rows ();
  int nc = columns ();

  switch (jtype)
    {
    case complex_scalar_constant:
    case scalar_constant:
      {
	int j = tree_to_mat_idx (tmp_j.double_value ());
	if (index_check (j, "column") < 0)
	  return tree_constant ();
	if (range_max_check (tree_to_mat_idx (ri.max ()), j, nr, nc) < 0)
	  return tree_constant ();
	retval = do_matrix_index (ri, j);
      }
      break;

    case complex_matrix_constant:
    case matrix_constant:
      {
	Matrix mj = tmp_j.matrix_value ();
	idx_vector jv (mj, user_pref.do_fortran_indexing, "column", nc);
	if (! jv)
	  return tree_constant ();

	if (jv.length () == 0)
	  {
	    Matrix mtmp;
	    retval = tree_constant (mtmp);
	  }
	else
	  {
	    if (range_max_check (tree_to_mat_idx (ri.max ()),
				 jv.max (), nr, nc) < 0)
	      return tree_constant ();
	    retval = do_matrix_index (ri, jv);
	  }
      }
      break;

    case string_constant:
      gripe_string_invalid ();
      break;

    case range_constant:
      {
	Range rj = tmp_j.range_value ();
	if (nc == 2 && is_zero_one (rj))
	  {
	    retval = do_matrix_index (ri, 1);
	  }
	else if (nc == 2 && is_one_zero (rj))
	  {
	    retval = do_matrix_index (ri, 0);
	  }
	else
	  {
	    if (index_check (rj, "column") < 0)
	      return tree_constant ();
	    if (range_max_check (tree_to_mat_idx (ri.max ()),
				 tree_to_mat_idx (rj.max ()), nr, nc) < 0)
	      return tree_constant ();
	    retval = do_matrix_index (ri, rj);
	  }
      }
      break;

    case magic_colon:
      {
	if (index_check (ri, "row") < 0)
	  return tree_constant ();
	if (range_max_check (tree_to_mat_idx (ri.max ()), 0, nr, nc) < 0)
	  return tree_constant ();
	retval = do_matrix_index (ri, magic_colon);
      }
      break;

    default:
      panic_impossible ();
      break;
    }

  return retval;
}

tree_constant
TC_REP::do_matrix_index (TC_REP::constant_type mci,
			 const tree_constant& j_arg) const
{
  tree_constant retval;

  tree_constant tmp_j = j_arg.make_numeric_or_range_or_magic ();

  if (error_state)
    return retval;

  TC_REP::constant_type jtype = tmp_j.const_type ();

  int nr = rows ();
  int nc = columns ();

  switch (jtype)
    {
    case complex_scalar_constant:
    case scalar_constant:
      {
	int j = tree_to_mat_idx (tmp_j.double_value ());
	if (j == -1 && nc == 1)
	  return Matrix ();
	if (index_check (j, "column") < 0)
	  return tree_constant ();
	if (range_max_check (0, j, nr, nc) < 0)
	  return tree_constant ();
	retval = do_matrix_index (magic_colon, j);
      }
      break;

    case complex_matrix_constant:
    case matrix_constant:
      {
	Matrix mj = tmp_j.matrix_value ();
	idx_vector jv (mj, user_pref.do_fortran_indexing, "column", nc);
	if (! jv)
	  return tree_constant ();

	if (jv.length () == 0)
	  {
	    Matrix mtmp;
	    retval = tree_constant (mtmp);
	  }
	else
	  {
	    if (range_max_check (0, jv.max (), nr, nc) < 0)
	      return tree_constant ();
	    retval = do_matrix_index (magic_colon, jv);
	  }
      }
      break;

    case string_constant:
      gripe_string_invalid ();
      break;

    case range_constant:
      {
	Range rj = tmp_j.range_value ();
	if (nc == 2 && is_zero_one (rj))
	  {
	    retval = do_matrix_index (magic_colon, 1);
	  }
	else if (nc == 2 && is_one_zero (rj))
	  {
	    retval = do_matrix_index (magic_colon, 0);
	  }
	else
	  {
	    if (index_check (rj, "column") < 0)
	      return tree_constant ();
	    if (range_max_check (0, tree_to_mat_idx (rj.max ()), nr, nc) < 0)
	      return tree_constant ();
	    retval = do_matrix_index (magic_colon, rj);
	  }
      }
      break;

    case magic_colon:
      retval = do_matrix_index (magic_colon, magic_colon);
      break;

    default:
      panic_impossible ();
      break;
    }

  return retval;
}

tree_constant
TC_REP::do_matrix_index (int i, int j) const
{
  tree_constant retval;

  if (type_tag == matrix_constant)
    retval = tree_constant (matrix->elem (i, j));
  else
    retval = tree_constant (complex_matrix->elem (i, j));

  return retval;
}

tree_constant
TC_REP::do_matrix_index (int i, const idx_vector& jv) const
{
  tree_constant retval;

  int jlen = jv.capacity ();

  CRMATRIX (m, cm, 1, jlen);

  for (int j = 0; j < jlen; j++)
    {
      int col = jv.elem (j);
      CRMATRIX_ASSIGN_REP_ELEM (m, cm, 0, j, i, col);
    }
  ASSIGN_CRMATRIX_TO (retval, m, cm);

  return retval;
}

tree_constant
TC_REP::do_matrix_index (int i, const Range& rj) const
{
  tree_constant retval;

  int jlen = rj.nelem ();

  CRMATRIX (m, cm, 1, jlen);

  double b = rj.base ();
  double increment = rj.inc ();
  for (int j = 0; j < jlen; j++)
    {
      double tmp = b + j * increment;
      int col = tree_to_mat_idx (tmp);
      CRMATRIX_ASSIGN_REP_ELEM (m, cm, 0, j, i, col);
    }

  ASSIGN_CRMATRIX_TO (retval, m, cm);

  return retval;
}

tree_constant
TC_REP::do_matrix_index (int i, TC_REP::constant_type mcj) const
{
  assert (mcj == magic_colon);

  tree_constant retval;

  int nc = columns ();

  CRMATRIX (m, cm, 1, nc);

  for (int j = 0; j < nc; j++)
    {
      CRMATRIX_ASSIGN_REP_ELEM (m, cm, 0, j, i, j);
    }

  ASSIGN_CRMATRIX_TO (retval, m, cm);

  return retval;
}

tree_constant
TC_REP::do_matrix_index (const idx_vector& iv, int j) const
{
  tree_constant retval;

  int ilen = iv.capacity ();

  CRMATRIX (m, cm, ilen, 1);

  for (int i = 0; i < ilen; i++)
    {
      int row = iv.elem (i);
      CRMATRIX_ASSIGN_REP_ELEM (m, cm, i, 0, row, j);
    }

  ASSIGN_CRMATRIX_TO (retval, m, cm);

  return retval;
}

tree_constant
TC_REP::do_matrix_index (const idx_vector& iv, const idx_vector& jv) const
{
  tree_constant retval;

  int ilen = iv.capacity ();
  int jlen = jv.capacity ();

  CRMATRIX (m, cm, ilen, jlen);

  for (int i = 0; i < ilen; i++)
    {
      int row = iv.elem (i);
      for (int j = 0; j < jlen; j++)
	{
	  int col = jv.elem (j);
	  CRMATRIX_ASSIGN_REP_ELEM (m, cm, i, j, row, col);
	}
    }

  ASSIGN_CRMATRIX_TO (retval, m, cm);

  return retval;
}

tree_constant
TC_REP::do_matrix_index (const idx_vector& iv, const Range& rj) const
{
  tree_constant retval;

  int ilen = iv.capacity ();
  int jlen = rj.nelem ();

  CRMATRIX (m, cm, ilen, jlen);

  double b = rj.base ();
  double increment = rj.inc ();

  for (int i = 0; i < ilen; i++)
    {
      int row = iv.elem (i);
      for (int j = 0; j < jlen; j++)
	{
	  double tmp = b + j * increment;
	  int col = tree_to_mat_idx (tmp);
	  CRMATRIX_ASSIGN_REP_ELEM (m, cm, i, j, row, col);
	}
    }

  ASSIGN_CRMATRIX_TO (retval, m, cm);

  return retval;
}

tree_constant
TC_REP::do_matrix_index (const idx_vector& iv,
			 TC_REP::constant_type mcj) const
{
  assert (mcj == magic_colon);

  tree_constant retval;

  int nc = columns ();
  int ilen = iv.capacity ();

  CRMATRIX (m, cm, ilen, nc);

  for (int j = 0; j < nc; j++)
    {
      for (int i = 0; i < ilen; i++)
	{
	  int row = iv.elem (i);
	  CRMATRIX_ASSIGN_REP_ELEM (m, cm, i, j, row, j);
	}
    }

  ASSIGN_CRMATRIX_TO (retval, m, cm);

  return retval;
}

tree_constant
TC_REP::do_matrix_index (const Range& ri, int j) const
{
  tree_constant retval;

  int ilen = ri.nelem ();

  CRMATRIX (m, cm, ilen, 1);

  double b = ri.base ();
  double increment = ri.inc ();
  for (int i = 0; i < ilen; i++)
    {
      double tmp = b + i * increment;
      int row = tree_to_mat_idx (tmp);
      CRMATRIX_ASSIGN_REP_ELEM (m, cm, i, 0, row, j);
    }

  ASSIGN_CRMATRIX_TO (retval, m, cm);

  return retval;
}

tree_constant
TC_REP::do_matrix_index (const Range& ri,
			 const idx_vector& jv) const
{
  tree_constant retval;

  int ilen = ri.nelem ();
  int jlen = jv.capacity ();

  CRMATRIX (m, cm, ilen, jlen);

  double b = ri.base ();
  double increment = ri.inc ();
  for (int i = 0; i < ilen; i++)
    {
      double tmp = b + i * increment;
      int row = tree_to_mat_idx (tmp);
      for (int j = 0; j < jlen; j++)
	{
	  int col = jv.elem (j);
	  CRMATRIX_ASSIGN_REP_ELEM (m, cm, i, j, row, col);
	}
    }

  ASSIGN_CRMATRIX_TO (retval, m, cm);

  return retval;
}

tree_constant
TC_REP::do_matrix_index (const Range& ri, const Range& rj) const
{
  tree_constant retval;

  int ilen = ri.nelem ();
  int jlen = rj.nelem ();

  CRMATRIX (m, cm, ilen, jlen);

  double ib = ri.base ();
  double iinc = ri.inc ();
  double jb = rj.base ();
  double jinc = rj.inc ();

  for (int i = 0; i < ilen; i++)
    {
      double itmp = ib + i * iinc;
      int row = tree_to_mat_idx (itmp);
      for (int j = 0; j < jlen; j++)
	{
	  double jtmp = jb + j * jinc;
	  int col = tree_to_mat_idx (jtmp);

	  CRMATRIX_ASSIGN_REP_ELEM (m, cm, i, j, row, col);
	}
    }

  ASSIGN_CRMATRIX_TO (retval, m, cm);

  return retval;
}

tree_constant
TC_REP::do_matrix_index (const Range& ri, TC_REP::constant_type mcj) const
{
  assert (mcj == magic_colon);

  tree_constant retval;

  int nc = columns ();

  int ilen = ri.nelem ();

  CRMATRIX (m, cm, ilen, nc);

  double ib = ri.base ();
  double iinc = ri.inc ();

  for (int i = 0; i < ilen; i++)
    {
      double itmp = ib + i * iinc;
      int row = tree_to_mat_idx (itmp);
      for (int j = 0; j < nc; j++)
	{
	  CRMATRIX_ASSIGN_REP_ELEM (m, cm, i, j, row, j);
	}
    }

  ASSIGN_CRMATRIX_TO (retval, m, cm);

  return retval;
}

tree_constant
TC_REP::do_matrix_index (TC_REP::constant_type mci, int j) const
{
  assert (mci == magic_colon);

  tree_constant retval;

  int nr = rows ();

  CRMATRIX (m, cm, nr, 1);

  for (int i = 0; i < nr; i++)
    {
      CRMATRIX_ASSIGN_REP_ELEM (m, cm, i, 0, i, j);
    }

  ASSIGN_CRMATRIX_TO (retval, m, cm);

  return retval;
}

tree_constant
TC_REP::do_matrix_index (TC_REP::constant_type mci,
			 const idx_vector& jv) const
{
  assert (mci == magic_colon);

  tree_constant retval;

  int nr = rows ();
  int jlen = jv.capacity ();

  CRMATRIX (m, cm, nr, jlen);

  for (int i = 0; i < nr; i++)
    {
      for (int j = 0; j < jlen; j++)
	{
	  int col = jv.elem (j);
	  CRMATRIX_ASSIGN_REP_ELEM (m, cm, i, j, i, col);
	}
    }

  ASSIGN_CRMATRIX_TO (retval, m, cm);

  return retval;
}

tree_constant
TC_REP::do_matrix_index (TC_REP::constant_type mci, const Range& rj) const
{
  assert (mci == magic_colon);

  tree_constant retval;

  int nr = rows ();
  int jlen = rj.nelem ();

  CRMATRIX (m, cm, nr, jlen);

  double jb = rj.base ();
  double jinc = rj.inc ();

  for (int j = 0; j < jlen; j++)
    {
      double jtmp = jb + j * jinc;
      int col = tree_to_mat_idx (jtmp);
      for (int i = 0; i < nr; i++)
	{
	  CRMATRIX_ASSIGN_REP_ELEM (m, cm, i, j, i, col);
	}
    }

  ASSIGN_CRMATRIX_TO (retval, m, cm);

  return retval;
}

tree_constant
TC_REP::do_matrix_index (TC_REP::constant_type mci,
			 TC_REP::constant_type mcj) const
{
  tree_constant retval;

  assert (mci == magic_colon && mcj == magic_colon);

  switch (type_tag)
    {
    case complex_scalar_constant:
      retval = *complex_scalar;
      break;

    case scalar_constant:
      retval = scalar;
      break;

    case complex_matrix_constant:
      retval = *complex_matrix;
      break;

    case matrix_constant:
      retval = *matrix;
      break;

    case range_constant:
      retval = *range;
      break;

    case string_constant:
      retval = string;
      break;

    case magic_colon:
    default:
      panic_impossible ();
      break;
    }

  return retval;
}

// -------------------------------------------------------------------
//
// Assignment operations for the tree-constant representation class.
//
// Leave the commented #includes below to make it easy to split this
// out again, should we want to do that.
//
// -------------------------------------------------------------------

// #ifdef HAVE_CONFIG_H
// #include <config.h>
// #endif

// #include <ctype.h>
// #include <string.h>
// #include <fstream.h>
// #include <iostream.h>
// #include <strstream.h>

// #include "mx-base.h"
// #include "Range.h"

// #include "arith-ops.h"
// #include "variables.h"
// #include "sysdep.h"
// #include "error.h"
// #include "gripes.h"
// #include "user-prefs.h"
// #include "utils.h"
// #include "pager.h"
// #include "pr-output.h"
// #include "tree-const.h"
// #include "idx-vector.h"
// #include "oct-map.h"

// #include "tc-inlines.h"

// Top-level tree-constant function that handles assignments.  Only
// decide if the left-hand side is currently a scalar or a matrix and
// hand off to other functions to do the real work.

void
TC_REP::assign (tree_constant& rhs, const Octave_object& args)
{
  tree_constant rhs_tmp = rhs.make_numeric ();

  if (error_state)
    return;

// This is easier than actually handling assignments to strings.
// An assignment to a range will normally require a conversion to a
// vector since it will normally destroy the equally-spaced property
// of the range elements.

  if (is_defined () && ! is_numeric_type ())
    force_numeric ();

  if (error_state)
    return;

  switch (type_tag)
    {
    case complex_scalar_constant:
    case scalar_constant:
    case unknown_constant:
      do_scalar_assignment (rhs_tmp, args);
      break;

    case complex_matrix_constant:
    case matrix_constant:
      do_matrix_assignment (rhs_tmp, args);
      break;

    default:
      ::error ("invalid assignment to %s", type_as_string ());
      break;
    }
}

// Assignments to scalars.  If resize_on_range_error is true,
// this can convert the left-hand side to a matrix.

void
TC_REP::do_scalar_assignment (const tree_constant& rhs,
			      const Octave_object& args) 
{
  assert (type_tag == unknown_constant
	  || type_tag == scalar_constant
	  || type_tag == complex_scalar_constant);

  int nargin = args.length ();

  if (rhs.is_zero_by_zero ())
    {
      if (valid_scalar_indices (args))
	{
	  if (type_tag == complex_scalar_constant)
	    delete complex_scalar;

	  matrix = new Matrix (0, 0);
	  type_tag = matrix_constant;
	}
      else if (! valid_zero_index (args))
	{
	  ::error ("invalid assigment of empty matrix to scalar");
	  return;
	}
    }
  else if (rhs.is_scalar_type () && valid_scalar_indices (args))
    {
      if (type_tag == unknown_constant || type_tag == scalar_constant)
	{
	  if (rhs.const_type () == scalar_constant)
	    {
	      scalar = rhs.double_value ();
	      type_tag = scalar_constant;
	    }
	  else if (rhs.const_type () == complex_scalar_constant)
	    {
	      complex_scalar = new Complex (rhs.complex_value ());
	      type_tag = complex_scalar_constant;
	    }
	  else
	    {
	      ::error ("invalid assignment to scalar");
	      return;
	    }
	}
      else
	{
	  if (rhs.const_type () == scalar_constant)
	    {
	      delete complex_scalar;
	      scalar = rhs.double_value ();
	      type_tag = scalar_constant;
	    }
	  else if (rhs.const_type () == complex_scalar_constant)
	    {
	      *complex_scalar = rhs.complex_value ();
	      type_tag = complex_scalar_constant;
	    }
	  else
	    {
	      ::error ("invalid assignment to scalar");
	      return;
	    }
	}
    }
  else if (user_pref.resize_on_range_error)
    {
      TC_REP::constant_type old_type_tag = type_tag;

      if (type_tag == complex_scalar_constant)
	{
	  Complex *old_complex = complex_scalar;
	  complex_matrix = new ComplexMatrix (1, 1, *complex_scalar);
	  type_tag = complex_matrix_constant;
	  delete old_complex;
	}
      else if (type_tag == scalar_constant)
	{
	  matrix = new Matrix (1, 1, scalar);
	  type_tag = matrix_constant;
	}

// If there is an error, the call to do_matrix_assignment should not
// destroy the current value.
// TC_REP::eval(int) will take
// care of converting single element matrices back to scalars.

      do_matrix_assignment (rhs, args);

// I don't think there's any other way to revert back to unknown
// constant types, so here it is.

      if (old_type_tag == unknown_constant && error_state)
	{
	  if (type_tag == matrix_constant)
	    delete matrix;
	  else if (type_tag == complex_matrix_constant)
	    delete complex_matrix;

	  type_tag = unknown_constant;
	}
    }
  else if (nargin > 2 || nargin < 1)
    ::error ("invalid index expression for scalar type");
  else
    ::error ("index invalid or out of range for scalar type");
}

// Assignments to matrices (and vectors).
//
// For compatibility with Matlab, we allow assignment of an empty
// matrix to an expression with empty indices to do nothing.

void
TC_REP::do_matrix_assignment (const tree_constant& rhs,
			      const Octave_object& args)
{
  assert (type_tag == unknown_constant
	  || type_tag == matrix_constant
	  || type_tag == complex_matrix_constant);

  if (type_tag == matrix_constant && rhs.is_complex_type ())
    {
      Matrix *old_matrix = matrix;
      complex_matrix = new ComplexMatrix (*matrix);
      type_tag = complex_matrix_constant;
      delete old_matrix;
    }
  else if (type_tag == unknown_constant)
    {
      if (rhs.is_complex_type ())
	{
	  complex_matrix = new ComplexMatrix ();
	  type_tag = complex_matrix_constant;
	}
      else
	{
	  matrix = new Matrix ();
	  type_tag = matrix_constant;
	}
    }

  int nargin = args.length ();

// The do_matrix_assignment functions can't handle empty matrices, so
// don't let any pass through here.
  switch (nargin)
    {
    case 1:
      {
	tree_constant arg = args(0);

	if (arg.is_undefined ())
	  ::error ("matrix index is undefined");
	else
	  do_matrix_assignment (rhs, arg);
      }
      break;

    case 2:
      {
	tree_constant arg_a = args(0);
	tree_constant arg_b = args(1);

	if (arg_a.is_undefined ())
	  ::error ("first matrix index is undefined");
	else if (arg_b.is_undefined ())
	  ::error ("second matrix index is undefined");
	else if (arg_a.is_empty () || arg_b.is_empty ())
	  {
	    if (! rhs.is_empty ())
	      {
		::error ("in assignment expression, a matrix index is empty");
		::error ("but the right hand side is not an empty matrix");
	      }
// XXX FIXME XXX -- to really be correct here, we should probably
// check to see if the assignment conforms, but that seems like more
// work than it's worth right now...
	  }
	else
	  do_matrix_assignment (rhs, arg_a, arg_b);
      }
      break;

    default:
      if (nargin == 0)
	::error ("matrix indices expected, but none provided");
      else
	::error ("too many indices for matrix expression");
      break;
    }
}

// Matrix assignments indexed by a single value.

void
TC_REP::do_matrix_assignment (const tree_constant& rhs,
			      const tree_constant& i_arg)
{
  int nr = rows ();
  int nc = columns ();

  if (user_pref.do_fortran_indexing || nr <= 1 || nc <= 1)
    {
      if (i_arg.is_empty ())
	{
	  if (! rhs.is_empty ())
	    {
	      ::error ("in assignment expression, matrix index is empty but");
	      ::error ("right hand side is not an empty matrix");
	    }
// XXX FIXME XXX -- to really be correct here, we should probably
// check to see if the assignment conforms, but that seems like more
// work than it's worth right now...

// The assignment functions can't handle empty matrices, so don't let
// any pass through here.
	  return;
	}

// We can't handle the case of assigning to a vector first, since even
// then, the two operations are not equivalent.  For example, the
// expression V(:) = M is handled differently depending on whether the
// user specified do_fortran_indexing = "true".

      if (user_pref.do_fortran_indexing)
	fortran_style_matrix_assignment (rhs, i_arg);
      else if (nr <= 1 || nc <= 1)
	vector_assignment (rhs, i_arg);
      else
	panic_impossible ();
    }
  else
    ::error ("single index only valid for row or column vector");
}

// Fortran-style assignments.  Matrices are assumed to be stored in
// column-major order and it is ok to use a single index for
// multi-dimensional matrices.

void
TC_REP::fortran_style_matrix_assignment (const tree_constant& rhs,
					 const tree_constant& i_arg)
{
  tree_constant tmp_i = i_arg.make_numeric_or_magic ();

  if (error_state)
    return;

  TC_REP::constant_type itype = tmp_i.const_type ();

  int nr = rows ();
  int nc = columns ();

  int rhs_nr = rhs.rows ();
  int rhs_nc = rhs.columns ();

  switch (itype)
    {
    case complex_scalar_constant:
    case scalar_constant:
      {
	double dval = tmp_i.double_value ();

	if (xisnan (dval))
	  {
	    error ("NaN is invalid as a matrix index");
	    return;
	  }

	int i = NINT (dval);
	int idx = i - 1;

	if (rhs_nr == 0 && rhs_nc == 0)
	  {
	    int len = nr * nc;

	    if (idx < len && len > 0)
	      {
		convert_to_row_or_column_vector ();

		nr = rows ();
		nc = columns ();

		if (nr == 1)
		  delete_column (idx);
		else if (nc == 1)
		  delete_row (idx);
		else
		  panic_impossible ();
	      }
	    else if (idx < 0)
	      {
		error ("invalid index = %d", idx+1);
	      }

	    return;
	  }

	if (index_check (idx, "") < 0)
	  return;

	if (nr <= 1 || nc <= 1)
	  {
	    maybe_resize (idx);
	    if (error_state)
	      return;
	  }
	else if (range_max_check (idx, nr * nc) < 0)
	  return;

	nr = rows ();
	nc = columns ();

	if (! indexed_assign_conforms (1, 1, rhs_nr, rhs_nc))
	  {
	    ::error ("for A(int) = X: X must be a scalar");
	    return;
	  }
	int ii = fortran_row (i, nr) - 1;
	int jj = fortran_column (i, nr) - 1;
	do_matrix_assignment (rhs, ii, jj);
      }
      break;

    case complex_matrix_constant:
    case matrix_constant:
      {
	Matrix mi = tmp_i.matrix_value ();
	int len = nr * nc;
	idx_vector ii (mi, 1, "", len);  // Always do fortran indexing here...
	if (! ii)
	  return;

	if (rhs_nr == 0 && rhs_nc == 0)
	  {
	    ii.sort_uniq ();
	    int num_to_delete = 0;
	    for (int i = 0; i < ii.length (); i++)
	      {
		if (ii.elem (i) < len)
		  num_to_delete++;
		else
		  break;
	      }

	    if (num_to_delete > 0)
	      {
		if (num_to_delete != ii.length ())
		  ii.shorten (num_to_delete);

		convert_to_row_or_column_vector ();

		nr = rows ();
		nc = columns ();

		if (nr == 1)
		  delete_columns (ii);
		else if (nc == 1)
		  delete_rows (ii);
		else
		  panic_impossible ();
	      }
	    return;
	  }

	if (nr <= 1 || nc <= 1)
	  {
	    maybe_resize (ii.max ());
	    if (error_state)
	      return;
	  }
	else if (range_max_check (ii.max (), len) < 0)
	  return;

	int ilen = ii.capacity ();

	if (ilen != rhs_nr * rhs_nc)
	  {
	    ::error ("A(matrix) = X: X and matrix must have the same number");
	    ::error ("of elements");
	  }
	else if (ilen == 1 && rhs.is_scalar_type ())
	  {
	    int nr = rows ();
	    int idx = ii.elem (0);
	    int ii = fortran_row (idx + 1, nr) - 1;
	    int jj = fortran_column (idx + 1, nr) - 1;

	    if (rhs.const_type () == scalar_constant)
	      matrix->elem (ii, jj) = rhs.double_value ();
	    else if (rhs.const_type () == complex_scalar_constant)
	      complex_matrix->elem (ii, jj) = rhs.complex_value ();
	    else
	      panic_impossible ();
	  }
	else
	  fortran_style_matrix_assignment (rhs, ii);
      }
      break;

    case string_constant:
      gripe_string_invalid ();
      break;

    case range_constant:
      gripe_range_invalid ();
      break;

    case magic_colon:
// a(:) = [] is equivalent to a(:,:) = [].
      if (rhs_nr == 0 && rhs_nc == 0)
	do_matrix_assignment (rhs, magic_colon, magic_colon);
      else
	fortran_style_matrix_assignment (rhs, magic_colon);
      break;

    default:
      panic_impossible ();
      break;
    }
}

// Fortran-style assignment for vector index.

void
TC_REP::fortran_style_matrix_assignment (const tree_constant& rhs,
					 idx_vector& i)
{
  assert (rhs.is_matrix_type ());

  int ilen = i.capacity ();

  REP_RHS_MATRIX (rhs, rhs_m, rhs_cm, rhs_nr, rhs_nc);

  int len = rhs_nr * rhs_nc;

  if (len == ilen)
    {
      int nr = rows ();
      if (rhs.const_type () == matrix_constant)
	{
	  double *cop_out = rhs_m.fortran_vec ();
	  for (int k = 0; k < len; k++)
	    {
	      int ii = fortran_row (i.elem (k) + 1, nr) - 1;
	      int jj = fortran_column (i.elem (k) + 1, nr) - 1;

	      matrix->elem (ii, jj) = *cop_out++;
	    }
	}
      else
	{
	  Complex *cop_out = rhs_cm.fortran_vec ();
	  for (int k = 0; k < len; k++)
	    {
	      int ii = fortran_row (i.elem (k) + 1, nr) - 1;
	      int jj = fortran_column (i.elem (k) + 1, nr) - 1;

	      complex_matrix->elem (ii, jj) = *cop_out++;
	    }
	}
    }
  else
    ::error ("number of rows and columns must match for indexed assignment");
}

// Fortran-style assignment for colon index.

void
TC_REP::fortran_style_matrix_assignment (const tree_constant& rhs,
					 TC_REP::constant_type mci)
{
  assert (rhs.is_matrix_type () && mci == TC_REP::magic_colon);

  int nr = rows ();
  int nc = columns ();

  REP_RHS_MATRIX (rhs, rhs_m, rhs_cm, rhs_nr, rhs_nc);

  int rhs_size = rhs_nr * rhs_nc;
  if (rhs_size == 0)
    {
      if (rhs.const_type () == matrix_constant)
	{
	  delete matrix;
	  matrix = new Matrix (0, 0);
	  return;
	}
      else
	panic_impossible ();
    }
  else if (nr*nc != rhs_size)
    {
      ::error ("A(:) = X: X and A must have the same number of elements");
      return;
    }

  if (rhs.const_type () == matrix_constant)
    {
      double *cop_out = rhs_m.fortran_vec ();
      for (int j = 0; j < nc; j++)
	for (int i = 0; i < nr; i++)
	  matrix->elem (i, j) = *cop_out++;
    }
  else
    {
      Complex *cop_out = rhs_cm.fortran_vec ();
      for (int j = 0; j < nc; j++)
	for (int i = 0; i < nr; i++)
	  complex_matrix->elem (i, j) = *cop_out++;
    }
}

// Assignments to vectors.  Hand off to other functions once we know
// what kind of index we have.  For a colon, it is the same as
// assignment to a matrix indexed by two colons.

void
TC_REP::vector_assignment (const tree_constant& rhs,
			   const tree_constant& i_arg)
{
  int nr = rows ();
  int nc = columns ();

  assert ((nr <= 1 || nc <= 1) && ! user_pref.do_fortran_indexing);

  tree_constant tmp_i = i_arg.make_numeric_or_range_or_magic ();

  if (error_state)
    return;

  TC_REP::constant_type itype = tmp_i.const_type ();

  switch (itype)
    {
    case complex_scalar_constant:
    case scalar_constant:
      {
	int i = tree_to_mat_idx (tmp_i.double_value ());
	if (index_check (i, "") < 0)
	  return;
	do_vector_assign (rhs, i);
      }
      break;

    case complex_matrix_constant:
    case matrix_constant:
      {
	Matrix mi = tmp_i.matrix_value ();
	int len = nr * nc;
	idx_vector iv (mi, user_pref.do_fortran_indexing, "", len);
	if (! iv)
	  return;

	do_vector_assign (rhs, iv);
      }
      break;

    case string_constant:
      gripe_string_invalid ();
      break;

    case range_constant:
      {
	Range ri = tmp_i.range_value ();
	int len = nr * nc;
	if (len == 2 && is_zero_one (ri))
	  {
	    do_vector_assign (rhs, 1);
	  }
	else if (len == 2 && is_one_zero (ri))
	  {
	    do_vector_assign (rhs, 0);
	  }
	else
	  {
	    if (index_check (ri, "") < 0)
	      return;
	    do_vector_assign (rhs, ri);
	  }
      }
      break;

    case magic_colon:
      {
	int rhs_nr = rhs.rows ();
	int rhs_nc = rhs.columns ();

	if (! indexed_assign_conforms (nr, nc, rhs_nr, rhs_nc))
	  {
	    ::error ("A(:) = X: X and A must have the same dimensions");
	    return;
	  }
	do_matrix_assignment (rhs, magic_colon, magic_colon);
      }
      break;

    default:
      panic_impossible ();
      break;
    }
}

// Check whether an indexed assignment to a vector is valid.

void
TC_REP::check_vector_assign (int rhs_nr, int rhs_nc, int ilen, const char *rm)
{
  int nr = rows ();
  int nc = columns ();

  if ((nr == 1 && nc == 1) || nr == 0 || nc == 0)  // No orientation.
    {
      if (! (ilen == rhs_nr || ilen == rhs_nc))
	{
	  ::error ("A(%s) = X: X and %s must have the same number of elements",
		 rm, rm);
	}
    }
  else if (nr == 1)  // Preserve current row orientation.
    {
      if (! (rhs_nr == 1 && rhs_nc == ilen))
	{
	  ::error ("A(%s) = X: where A is a row vector, X must also be a", rm);
	  ::error ("row vector with the same number of elements as %s", rm);
	}
    }
  else if (nc == 1)  // Preserve current column orientation.
    {
      if (! (rhs_nc == 1 && rhs_nr == ilen))
	{
	  ::error ("A(%s) = X: where A is a column vector, X must also be", rm);
	  ::error ("a column vector with the same number of elements as %s", rm);
	}
    }
  else
    panic_impossible ();
}

// Assignment to a vector with an integer index.

void
TC_REP::do_vector_assign (const tree_constant& rhs, int i)
{
  int rhs_nr = rhs.rows ();
  int rhs_nc = rhs.columns ();

  if (indexed_assign_conforms (1, 1, rhs_nr, rhs_nc))
    {
      maybe_resize (i);
      if (error_state)
	return;

      int nr = rows ();
      int nc = columns ();

      if (nr == 1)
	{
	  REP_ELEM_ASSIGN (0, i, rhs.double_value (), rhs.complex_value (),
			   rhs.is_real_type ());
	}
      else if (nc == 1)
	{
	  REP_ELEM_ASSIGN (i, 0, rhs.double_value (), rhs.complex_value (),
			   rhs.is_real_type ());
	}
      else
	panic_impossible ();
    }
  else if (rhs_nr == 0 && rhs_nc == 0)
    {
      int nr = rows ();
      int nc = columns ();

      int len = MAX (nr, nc);

      if (i < 0 || i >= len || (nr == 0 && nc == 0))
	{
	  ::error ("A(int) = []: index out of range");
	  return;
	}

      if (nr == 0 && nc > 0)
	resize (0, nc - 1);
      else if (nc == 0 && nr > 0)
	resize (nr - 1, 0);
      else if (nr == 1)
	delete_column (i);
      else if (nc == 1)
	delete_row (i);
      else
	panic_impossible ();
    }
  else
    {
      ::error ("for A(int) = X: X must be a scalar");
      return;
    }
}

// Assignment to a vector with a vector index.

void
TC_REP::do_vector_assign (const tree_constant& rhs, idx_vector& iv)
{
  if (rhs.is_zero_by_zero ())
    {
      int nr = rows ();
      int nc = columns ();

      int len = MAX (nr, nc);

      if (iv.max () >= len)
	{
	  ::error ("A(matrix) = []: index out of range");
	  return;
	}

      if (nr == 1)
	delete_columns (iv);
      else if (nc == 1)
	delete_rows (iv);
      else
	panic_impossible ();
    }
  else if (rhs.is_scalar_type ())
    {
      int nr = rows ();
      int nc = columns ();

      if (iv.capacity () == 1)
	{
	  int idx = iv.elem (0);

	  if (nr == 1)
	    {
	      REP_ELEM_ASSIGN (0, idx, rhs.double_value (),
			       rhs.complex_value (), rhs.is_real_type ());
	    }
	  else if (nc == 1)
	    {
	      REP_ELEM_ASSIGN (idx, 0, rhs.double_value (),
			       rhs.complex_value (), rhs.is_real_type ());
	    }
	  else
	    panic_impossible ();
	}
      else
	{
	  if (nr == 1)
	    {
	      ::error ("A(matrix) = X: where A is a row vector, X must also be a");
	      ::error ("row vector with the same number of elements as matrix");
	    }
	  else if (nc == 1)
	    {
	      ::error ("A(matrix) = X: where A is a column vector, X must also be a");
	      ::error ("column vector with the same number of elements as matrix");
	    }
	  else
	    panic_impossible ();
	}
    }
  else if (rhs.is_matrix_type ())
    {
      REP_RHS_MATRIX (rhs, rhs_m, rhs_cm, rhs_nr, rhs_nc);

      int ilen = iv.capacity ();
      check_vector_assign (rhs_nr, rhs_nc, ilen, "matrix");
      if (error_state)
	return;

      force_orient f_orient = no_orient;
      if (rhs_nr == 1 && rhs_nc != 1)
	f_orient = row_orient;
      else if (rhs_nc == 1 && rhs_nr != 1)
	f_orient = column_orient;

      maybe_resize (iv.max (), f_orient);
      if (error_state)
	return;

      int nr = rows ();
      int nc = columns ();

      if (nr == 1 && rhs_nr == 1)
	{
	  for (int i = 0; i < iv.capacity (); i++)
	    REP_ELEM_ASSIGN (0, iv.elem (i), rhs_m.elem (0, i),
			     rhs_cm.elem (0, i), rhs.is_real_type ());
	}
      else if (nc == 1 && rhs_nc == 1)
	{
	  for (int i = 0; i < iv.capacity (); i++)
	    REP_ELEM_ASSIGN (iv.elem (i), 0, rhs_m.elem (i, 0),
			     rhs_cm.elem (i, 0), rhs.is_real_type ());
	}
      else
	::error ("A(vector) = X: X must be the same size as vector");
    }
  else
    panic_impossible ();
}

// Assignment to a vector with a range index.

void
TC_REP::do_vector_assign (const tree_constant& rhs, Range& ri)
{
  if (rhs.is_zero_by_zero ())
    {
      int nr = rows ();
      int nc = columns ();

      int len = MAX (nr, nc);

      int b = tree_to_mat_idx (ri.min ());
      int l = tree_to_mat_idx (ri.max ());
      if (b < 0 || l >= len)
	{
	  ::error ("A(range) = []: index out of range");
	  return;
	}

      if (nr == 1)
	delete_columns (ri);
      else if (nc == 1)
	delete_rows (ri);
      else
	panic_impossible ();
    }
  else if (rhs.is_scalar_type ())
    {
      int nr = rows ();
      int nc = columns ();

      if (nr == 1)
	{
	  ::error ("A(range) = X: where A is a row vector, X must also be a");
	  ::error ("row vector with the same number of elements as range");
	}
      else if (nc == 1)
	{
	  ::error ("A(range) = X: where A is a column vector, X must also be a");
	  ::error ("column vector with the same number of elements as range");
	}
      else
	panic_impossible ();
    }
  else if (rhs.is_matrix_type ())
    {
      REP_RHS_MATRIX (rhs, rhs_m, rhs_cm, rhs_nr, rhs_nc);

      int ilen = ri.nelem ();
      check_vector_assign (rhs_nr, rhs_nc, ilen, "range");
      if (error_state)
	return;

      force_orient f_orient = no_orient;
      if (rhs_nr == 1 && rhs_nc != 1)
	f_orient = row_orient;
      else if (rhs_nc == 1 && rhs_nr != 1)
	f_orient = column_orient;

      maybe_resize (tree_to_mat_idx (ri.max ()), f_orient);
      if (error_state)
	return;

      int nr = rows ();
      int nc = columns ();

      double b = ri.base ();
      double increment = ri.inc ();

      if (nr == 1)
	{
	  for (int i = 0; i < ri.nelem (); i++)
	    {
	      double tmp = b + i * increment;
	      int col = tree_to_mat_idx (tmp);
	      REP_ELEM_ASSIGN (0, col, rhs_m.elem (0, i), rhs_cm.elem (0, i),
			       rhs.is_real_type ());
	    }
	}
      else if (nc == 1)
	{
	  for (int i = 0; i < ri.nelem (); i++)
	    {
	      double tmp = b + i * increment;
	      int row = tree_to_mat_idx (tmp);
	      REP_ELEM_ASSIGN (row, 0, rhs_m.elem (i, 0), rhs_cm.elem (i, 0),
			       rhs.is_real_type ());
	    }
	}
      else
	panic_impossible ();
    }
  else
    panic_impossible ();
}

// Matrix assignment indexed by two values.  This function determines
// the type of the first arugment, checks as much as possible, and
// then calls one of a set of functions to handle the specific cases:
//
//   M (integer, arg2) = RHS  (MA1)
//   M (vector,  arg2) = RHS  (MA2)
//   M (range,   arg2) = RHS  (MA3)
//   M (colon,   arg2) = RHS  (MA4)
//
// Each of those functions determines the type of the second argument
// and calls another function to handle the real work of doing the
// assignment.

void
TC_REP::do_matrix_assignment (const tree_constant& rhs,
			      const tree_constant& i_arg,
			      const tree_constant& j_arg)
{
  tree_constant tmp_i = i_arg.make_numeric_or_range_or_magic ();

  if (error_state)
    return;

  TC_REP::constant_type itype = tmp_i.const_type ();

  switch (itype)
    {
    case complex_scalar_constant:
    case scalar_constant:
      {
	int i = tree_to_mat_idx (tmp_i.double_value ());
	do_matrix_assignment (rhs, i, j_arg);
      }
      break;

    case complex_matrix_constant:
    case matrix_constant:
      {
	Matrix mi = tmp_i.matrix_value ();
	idx_vector iv (mi, user_pref.do_fortran_indexing, "row", rows ());
	if (! iv)
	  return;

	do_matrix_assignment (rhs, iv, j_arg);
      }
      break;

    case string_constant:
      gripe_string_invalid ();
      break;

    case range_constant:
      {
	Range ri = tmp_i.range_value ();
	int nr = rows ();
	if (nr == 2 && is_zero_one (ri))
	  {
	    do_matrix_assignment (rhs, 1, j_arg);
	  }
	else if (nr == 2 && is_one_zero (ri))
	  {
	    do_matrix_assignment (rhs, 0, j_arg);
	  }
	else
	  {
	    if (index_check (ri, "row") < 0)
	      return;
	    do_matrix_assignment (rhs, ri, j_arg);
	  }
      }
      break;

    case magic_colon:
      do_matrix_assignment (rhs, magic_colon, j_arg);
      break;

    default:
      panic_impossible ();
      break;
    }
}

/* MA1 */
void
TC_REP::do_matrix_assignment (const tree_constant& rhs, int i,
			      const tree_constant& j_arg)
{
  tree_constant tmp_j = j_arg.make_numeric_or_range_or_magic ();

  if (error_state)
    return;

  TC_REP::constant_type jtype = tmp_j.const_type ();

  int rhs_nr = rhs.rows ();
  int rhs_nc = rhs.columns ();

  switch (jtype)
    {
    case complex_scalar_constant:
    case scalar_constant:
      {
	if (index_check (i, "row") < 0)
	  return;
	int j = tree_to_mat_idx (tmp_j.double_value ());
	if (index_check (j, "column") < 0)
	  return;
	if (! indexed_assign_conforms (1, 1, rhs_nr, rhs_nc))
	  {
	    ::error ("A(int,int) = X, X must be a scalar");
	    return;
	  }
	maybe_resize (i, j);
	if (error_state)
	  return;

	do_matrix_assignment (rhs, i, j);
      }
      break;

    case complex_matrix_constant:
    case matrix_constant:
      {
	if (index_check (i, "row") < 0)
	  return;
	Matrix mj = tmp_j.matrix_value ();
	idx_vector jv (mj, user_pref.do_fortran_indexing, "column",
		       columns ());
	if (! jv)
	  return;

	if (! indexed_assign_conforms (1, jv.capacity (), rhs_nr, rhs_nc))
	  {
	    ::error ("A(int,matrix) = X: X must be a row vector with the same");
	    ::error ("number of elements as matrix");
	    return;
	  }
	maybe_resize (i, jv.max ());
	if (error_state)
	  return;

	do_matrix_assignment (rhs, i, jv);
      }
      break;

    case string_constant:
      gripe_string_invalid ();
      break;

    case range_constant:
      {
	if (index_check (i, "row") < 0)
	  return;
	Range rj = tmp_j.range_value ();
	if (! indexed_assign_conforms (1, rj.nelem (), rhs_nr, rhs_nc))
	  {
	    ::error ("A(int,range) = X: X must be a row vector with the same");
	    ::error ("number of elements as range");
	    return;
	  }

	int nc = columns ();
	if (nc == 2 && is_zero_one (rj) && rhs_nc == 1)
	  {
	    do_matrix_assignment (rhs, i, 1);
	  }
	else if (nc == 2 && is_one_zero (rj) && rhs_nc == 1)
	  {
	    do_matrix_assignment (rhs, i, 0);
	  }
	else
	  {
	    if (index_check (rj, "column") < 0)
	      return;
	    maybe_resize (i, tree_to_mat_idx (rj.max ()));
	    if (error_state)
	      return;

	    do_matrix_assignment (rhs, i, rj);
	  }
      }
      break;

    case magic_colon:
      {
	int nc = columns ();
	int nr = rows ();
	if (i == -1 && nr == 1 && rhs_nr == 0 && rhs_nc == 0
	    || index_check (i, "row") < 0)
	  return;
	else if (nc == 0 && nr == 0 && rhs_nr == 1)
	  {
	    if (rhs.is_complex_type ())
	      {
		complex_matrix = new ComplexMatrix ();
		type_tag = complex_matrix_constant;
	      }
	    else
	      {
		matrix = new Matrix ();
		type_tag = matrix_constant;
	      }
	    maybe_resize (i, rhs_nc-1);
	    if (error_state)
	      return;
	  }
	else if (indexed_assign_conforms (1, nc, rhs_nr, rhs_nc))
	  {
	    maybe_resize (i, nc-1);
	    if (error_state)
	      return;
	  }
	else if (rhs_nr == 0 && rhs_nc == 0)
	  {
	    if (i < 0 || i >= nr)
	      {
		::error ("A(int,:) = []: row index out of range");
		return;
	      }
	  }
	else
	  {
	    ::error ("A(int,:) = X: X must be a row vector with the same");
	    ::error ("number of columns as A");
	    return;
	  }

	do_matrix_assignment (rhs, i, magic_colon);
      }
      break;

    default:
      panic_impossible ();
      break;
    }
}

/* MA2 */
void
TC_REP::do_matrix_assignment (const tree_constant& rhs,
			      idx_vector& iv, const tree_constant& j_arg)
{
  tree_constant tmp_j = j_arg.make_numeric_or_range_or_magic ();

  if (error_state)
    return;

  TC_REP::constant_type jtype = tmp_j.const_type ();

  int rhs_nr = rhs.rows ();
  int rhs_nc = rhs.columns ();

  switch (jtype)
    {
    case complex_scalar_constant:
    case scalar_constant:
      {
	int j = tree_to_mat_idx (tmp_j.double_value ());
	if (index_check (j, "column") < 0)
	  return;
	if (! indexed_assign_conforms (iv.capacity (), 1, rhs_nr, rhs_nc))
	  {
	    ::error ("A(matrix,int) = X: X must be a column vector with the");
	    ::error ("same number of elements as matrix");
	    return;
	  }
	maybe_resize (iv.max (), j);
	if (error_state)
	  return;

	do_matrix_assignment (rhs, iv, j);
      }
      break;

    case complex_matrix_constant:
    case matrix_constant:
      {
	Matrix mj = tmp_j.matrix_value ();
	idx_vector jv (mj, user_pref.do_fortran_indexing, "column",
		       columns ());
	if (! jv)
	  return;

	if (! indexed_assign_conforms (iv.capacity (), jv.capacity (),
				       rhs_nr, rhs_nc))
	  {
	    ::error ("A(r_mat,c_mat) = X: the number of rows in X must match");
	    ::error ("the number of elements in r_mat and the number of");
	    ::error ("columns in X must match the number of elements in c_mat");
	    return;
	  }
	maybe_resize (iv.max (), jv.max ());
	if (error_state)
	  return;

	do_matrix_assignment (rhs, iv, jv);
      }
      break;

    case string_constant:
      gripe_string_invalid ();
      break;

    case range_constant:
      {
	Range rj = tmp_j.range_value ();
	if (! indexed_assign_conforms (iv.capacity (), rj.nelem (),
				       rhs_nr, rhs_nc))
	  {
	    ::error ("A(matrix,range) = X: the number of rows in X must match");
	    ::error ("the number of elements in matrix and the number of");
	    ::error ("columns in X must match the number of elements in range");
	    return;
	  }

	int nc = columns ();
	if (nc == 2 && is_zero_one (rj) && rhs_nc == 1)
	  {
	    do_matrix_assignment (rhs, iv, 1);
	  }
	else if (nc == 2 && is_one_zero (rj) && rhs_nc == 1)
	  {
	    do_matrix_assignment (rhs, iv, 0);
	  }
	else
	  {
	    if (index_check (rj, "column") < 0)
	      return;
	    maybe_resize (iv.max (), tree_to_mat_idx (rj.max ()));
	    if (error_state)
	      return;

	    do_matrix_assignment (rhs, iv, rj);
	  }
      }
      break;

    case magic_colon:
      {
	int nc = columns ();
	int new_nc = nc;
	if (nc == 0)
	  new_nc = rhs_nc;

	if (indexed_assign_conforms (iv.capacity (), new_nc,
				     rhs_nr, rhs_nc))
	  {
	    maybe_resize (iv.max (), new_nc-1);
	    if (error_state)
	      return;
	  }
	else if (rhs_nr == 0 && rhs_nc == 0)
	  {
	    if (iv.max () >= rows ())
	      {
		::error ("A(matrix,:) = []: row index out of range");
		return;
	      }
	  }
	else
	  {
	    ::error ("A(matrix,:) = X: the number of rows in X must match the");
	    ::error ("number of elements in matrix, and the number of columns");
	    ::error ("in X must match the number of columns in A");
	    return;
	  }

	do_matrix_assignment (rhs, iv, magic_colon);
      }
      break;

    default:
      panic_impossible ();
      break;
    }
}

/* MA3 */
void
TC_REP::do_matrix_assignment (const tree_constant& rhs, Range& ri,
			      const tree_constant& j_arg)
{
  tree_constant tmp_j = j_arg.make_numeric_or_range_or_magic ();

  if (error_state)
    return;

  TC_REP::constant_type jtype = tmp_j.const_type ();

  int rhs_nr = rhs.rows ();
  int rhs_nc = rhs.columns ();

  switch (jtype)
    {
    case complex_scalar_constant:
    case scalar_constant:
      {
	int j = tree_to_mat_idx (tmp_j.double_value ());
	if (index_check (j, "column") < 0)
	  return;
	if (! indexed_assign_conforms (ri.nelem (), 1, rhs_nr, rhs_nc))
	  {
	    ::error ("A(range,int) = X: X must be a column vector with the");
	    ::error ("same number of elements as range");
	    return;
	  }
	maybe_resize (tree_to_mat_idx (ri.max ()), j);
	if (error_state)
	  return;

	do_matrix_assignment (rhs, ri, j);
      }
      break;

    case complex_matrix_constant:
    case matrix_constant:
      {
	Matrix mj = tmp_j.matrix_value ();
	idx_vector jv (mj, user_pref.do_fortran_indexing, "column",
		       columns ());
	if (! jv)
	  return;

	if (! indexed_assign_conforms (ri.nelem (), jv.capacity (),
				       rhs_nr, rhs_nc))
	  {
	    ::error ("A(range,matrix) = X: the number of rows in X must match");
	    ::error ("the number of elements in range and the number of");
	    ::error ("columns in X must match the number of elements in matrix");
	    return;
	  }
	maybe_resize (tree_to_mat_idx (ri.max ()), jv.max ());
	if (error_state)
	  return;

	do_matrix_assignment (rhs, ri, jv);
      }
      break;

    case string_constant:
      gripe_string_invalid ();
      break;

    case range_constant:
      {
	Range rj = tmp_j.range_value ();
	if (! indexed_assign_conforms (ri.nelem (), rj.nelem (),
				       rhs_nr, rhs_nc))
	  {
	    ::error ("A(r_range,c_range) = X: the number of rows in X must");
	    ::error ("match the number of elements in r_range and the number");
	    ::error ("of columns in X must match the number of elements in");
	    ::error ("c_range");
	    return;
	  }

	int nc = columns ();
	if (nc == 2 && is_zero_one (rj) && rhs_nc == 1)
	  {
	    do_matrix_assignment (rhs, ri, 1);
	  }
	else if (nc == 2 && is_one_zero (rj) && rhs_nc == 1)
	  {
	    do_matrix_assignment (rhs, ri, 0);
	  }
	else
	  {
	    if (index_check (rj, "column") < 0)
	      return;

	    maybe_resize (tree_to_mat_idx (ri.max ()),
			  tree_to_mat_idx (rj.max ()));

	    if (error_state)
	      return;

	    do_matrix_assignment (rhs, ri, rj);
	  }
      }
      break;

    case magic_colon:
      {
	int nc = columns ();
	int new_nc = nc;
	if (nc == 0)
	  new_nc = rhs_nc;

	if (indexed_assign_conforms (ri.nelem (), new_nc, rhs_nr, rhs_nc))
	  {
	    maybe_resize (tree_to_mat_idx (ri.max ()), new_nc-1);
	    if (error_state)
	      return;
	  }
	else if (rhs_nr == 0 && rhs_nc == 0)
	  {
	    int b = tree_to_mat_idx (ri.min ());
	    int l = tree_to_mat_idx (ri.max ());
	    if (b < 0 || l >= rows ())
	      {
		::error ("A(range,:) = []: row index out of range");
		return;
	      }
	  }
	else
	  {
	    ::error ("A(range,:) = X: the number of rows in X must match the");
	    ::error ("number of elements in range, and the number of columns");
	    ::error ("in X must match the number of columns in A");
	    return;
	  }

	do_matrix_assignment (rhs, ri, magic_colon);
      }
      break;

    default:
      panic_impossible ();
      break;
    }
}

/* MA4 */
void
TC_REP::do_matrix_assignment (const tree_constant& rhs,
			      TC_REP::constant_type i,
			      const tree_constant& j_arg)
{
  tree_constant tmp_j = j_arg.make_numeric_or_range_or_magic ();

  if (error_state)
    return;

  TC_REP::constant_type jtype = tmp_j.const_type ();

  int rhs_nr = rhs.rows ();
  int rhs_nc = rhs.columns ();

  switch (jtype)
    {
    case complex_scalar_constant:
    case scalar_constant:
      {
	int j = tree_to_mat_idx (tmp_j.double_value ());
	int nr = rows ();
	int nc = columns ();
	if (j == -1 && nc == 1 && rhs_nr == 0 && rhs_nc == 0
	    || index_check (j, "column") < 0)
	  return;
	if (nr == 0 && nc == 0 && rhs_nc == 1)
	  {
	    if (rhs.is_complex_type ())
	      {
		complex_matrix = new ComplexMatrix ();
		type_tag = complex_matrix_constant;
	      }
	    else
	      {
		matrix = new Matrix ();
		type_tag = matrix_constant;
	      }
	    maybe_resize (rhs_nr-1, j);
	    if (error_state)
	      return;
	  }
	else if (indexed_assign_conforms (nr, 1, rhs_nr, rhs_nc))
	  {
	    maybe_resize (nr-1, j);
	    if (error_state)
	      return;
	  }
	else if (rhs_nr == 0 && rhs_nc == 0)
	  {
	    if (j < 0 || j >= nc)
	      {
		::error ("A(:,int) = []: column index out of range");
		return;
	      }
	  }
	else
	  {
	    ::error ("A(:,int) = X: X must be a column vector with the same");
	    ::error ("number of rows as A");
	    return;
	  }

	do_matrix_assignment (rhs, magic_colon, j);
      }
      break;

    case complex_matrix_constant:
    case matrix_constant:
      {
	Matrix mj = tmp_j.matrix_value ();
	idx_vector jv (mj, user_pref.do_fortran_indexing, "column",
		       columns ());
	if (! jv)
	  return;

	int nr = rows ();
	int new_nr = nr;
	if (nr == 0)
	  new_nr = rhs_nr;

	if (indexed_assign_conforms (new_nr, jv.capacity (),
				     rhs_nr, rhs_nc))
	  {
	    maybe_resize (new_nr-1, jv.max ());
	    if (error_state)
	      return;
	  }
	else if (rhs_nr == 0 && rhs_nc == 0)
	  {
	    if (jv.max () >= columns ())
	      {
		::error ("A(:,matrix) = []: column index out of range");
		return;
	      }
	  }
	else
	  {
	    ::error ("A(:,matrix) = X: the number of rows in X must match the");
	    ::error ("number of rows in A, and the number of columns in X must");
	    ::error ("match the number of elements in matrix");
	    return;
	  }

	do_matrix_assignment (rhs, magic_colon, jv);
      }
      break;

    case string_constant:
      gripe_string_invalid ();
      break;

    case range_constant:
      {
	Range rj = tmp_j.range_value ();
	int nr = rows ();
	int new_nr = nr;
	if (nr == 0)
	  new_nr = rhs_nr;

	if (indexed_assign_conforms (new_nr, rj.nelem (), rhs_nr, rhs_nc))
	  {
	    int nc = columns ();
	    if (nc == 2 && is_zero_one (rj) && rhs_nc == 1)
	      {
		do_matrix_assignment (rhs, magic_colon, 1);
	      }
	    else if (nc == 2 && is_one_zero (rj) && rhs_nc == 1)
	      {
		do_matrix_assignment (rhs, magic_colon, 0);
	      }
	    else
	      {
		if (index_check (rj, "column") < 0)
		  return;
		maybe_resize (new_nr-1, tree_to_mat_idx (rj.max ()));
		if (error_state)
		  return;
	      }
	  }
	else if (rhs_nr == 0 && rhs_nc == 0)
	  {
	    int b = tree_to_mat_idx (rj.min ());
	    int l = tree_to_mat_idx (rj.max ());
	    if (b < 0 || l >= columns ())
	      {
		::error ("A(:,range) = []: column index out of range");
		return;
	      }
	  }
	else
	  {
	    ::error ("A(:,range) = X: the number of rows in X must match the");
	    ::error ("number of rows in A, and the number of columns in X");
	    ::error ("must match the number of elements in range");
	    return;
	  }

	do_matrix_assignment (rhs, magic_colon, rj);
      }
      break;

    case magic_colon:
// a(:,:) = foo is equivalent to a = foo.
      do_matrix_assignment (rhs, magic_colon, magic_colon);
      break;

    default:
      panic_impossible ();
      break;
    }
}

// Functions that actually handle assignment to a matrix using two
// index values.
//
//                   idx2
//            +---+---+----+----+
//   idx1     | i | v |  r | c  |
//   ---------+---+---+----+----+
//   integer  | 1 | 5 |  9 | 13 |
//   ---------+---+---+----+----+
//   vector   | 2 | 6 | 10 | 14 |
//   ---------+---+---+----+----+
//   range    | 3 | 7 | 11 | 15 |
//   ---------+---+---+----+----+
//   colon    | 4 | 8 | 12 | 16 |
//   ---------+---+---+----+----+

/* 1 */
void
TC_REP::do_matrix_assignment (const tree_constant& rhs, int i, int j)
{
  REP_ELEM_ASSIGN (i, j, rhs.double_value (), rhs.complex_value (),
		   rhs.is_real_type ());
}

/* 2 */
void
TC_REP::do_matrix_assignment (const tree_constant& rhs, int i, idx_vector& jv)
{
  REP_RHS_MATRIX (rhs, rhs_m, rhs_cm, rhs_nr, rhs_nc);

  for (int j = 0; j < jv.capacity (); j++)
    REP_ELEM_ASSIGN (i, jv.elem (j), rhs_m.elem (0, j),
		     rhs_cm.elem (0, j), rhs.is_real_type ());
}

/* 3 */
void
TC_REP::do_matrix_assignment (const tree_constant& rhs, int i, Range& rj)
{
  REP_RHS_MATRIX (rhs, rhs_m, rhs_cm, rhs_nr, rhs_nc);

  double b = rj.base ();
  double increment = rj.inc ();

  for (int j = 0; j < rj.nelem (); j++)
    {
      double tmp = b + j * increment;
      int col = tree_to_mat_idx (tmp);
      REP_ELEM_ASSIGN (i, col, rhs_m.elem (0, j), rhs_cm.elem (0, j),
		       rhs.is_real_type ());
    }
}

/* 4 */
void
TC_REP::do_matrix_assignment (const tree_constant& rhs, int i,
			      TC_REP::constant_type mcj)
{
  assert (mcj == magic_colon);

  int nc = columns ();

  if (rhs.is_zero_by_zero ())
    {
      delete_row (i);
    }
  else if (rhs.is_matrix_type ())
    {
      REP_RHS_MATRIX (rhs, rhs_m, rhs_cm, rhs_nr, rhs_nc);

      for (int j = 0; j < nc; j++)
	REP_ELEM_ASSIGN (i, j, rhs_m.elem (0, j), rhs_cm.elem (0, j),
			 rhs.is_real_type ());
    }
  else if (rhs.is_scalar_type () && nc == 1)
    {
      REP_ELEM_ASSIGN (i, 0, rhs.double_value (),
		       rhs.complex_value (), rhs.is_real_type ());
    }
  else
    panic_impossible ();
}

/* 5 */
void
TC_REP::do_matrix_assignment (const tree_constant& rhs,
			      idx_vector& iv, int j)
{
  REP_RHS_MATRIX (rhs, rhs_m, rhs_cm, rhs_nr, rhs_nc);

  for (int i = 0; i < iv.capacity (); i++)
    {
      int row = iv.elem (i);
      REP_ELEM_ASSIGN (row, j, rhs_m.elem (i, 0),
		       rhs_cm.elem (i, 0), rhs.is_real_type ());
    }
}

/* 6 */
void
TC_REP::do_matrix_assignment (const tree_constant& rhs,
			      idx_vector& iv, idx_vector& jv)
{
  REP_RHS_MATRIX (rhs, rhs_m, rhs_cm, rhs_nr, rhs_nc);

  for (int i = 0; i < iv.capacity (); i++)
    {
      int row = iv.elem (i);
      for (int j = 0; j < jv.capacity (); j++)
	{
	  int col = jv.elem (j);
	  REP_ELEM_ASSIGN (row, col, rhs_m.elem (i, j),
			   rhs_cm.elem (i, j), rhs.is_real_type ());
	}
    }
}

/* 7 */
void
TC_REP::do_matrix_assignment (const tree_constant& rhs,
			      idx_vector& iv, Range& rj)
{
  REP_RHS_MATRIX (rhs, rhs_m, rhs_cm, rhs_nr, rhs_nc);

  double b = rj.base ();
  double increment = rj.inc ();

  for (int i = 0; i < iv.capacity (); i++)
    {
      int row = iv.elem (i);
      for (int j = 0; j < rj.nelem (); j++)
	{
	  double tmp = b + j * increment;
	  int col = tree_to_mat_idx (tmp);
	  REP_ELEM_ASSIGN (row, col, rhs_m.elem (i, j),
			   rhs_cm.elem (i, j), rhs.is_real_type ());
	}
    }
}

/* 8 */
void
TC_REP::do_matrix_assignment (const tree_constant& rhs,
			      idx_vector& iv, TC_REP::constant_type mcj)
{
  assert (mcj == magic_colon);

  if (rhs.is_zero_by_zero ())
    {
      delete_rows (iv);
    }
  else
    {
      REP_RHS_MATRIX (rhs, rhs_m, rhs_cm, rhs_nr, rhs_nc);

      int nc = columns ();

      for (int j = 0; j < nc; j++)
	{
	  for (int i = 0; i < iv.capacity (); i++)
	    {
	      int row = iv.elem (i);
	      REP_ELEM_ASSIGN (row, j, rhs_m.elem (i, j),
			       rhs_cm.elem (i, j), rhs.is_real_type ());
	    }
	}
    }
}

/* 9 */
void
TC_REP::do_matrix_assignment (const tree_constant& rhs, Range& ri, int j)
{
  REP_RHS_MATRIX (rhs, rhs_m, rhs_cm, rhs_nr, rhs_nc);

  double b = ri.base ();
  double increment = ri.inc ();

  for (int i = 0; i < ri.nelem (); i++)
    {
      double tmp = b + i * increment;
      int row = tree_to_mat_idx (tmp);
      REP_ELEM_ASSIGN (row, j, rhs_m.elem (i, 0),
		       rhs_cm.elem (i, 0), rhs.is_real_type ());
    }
}

/* 10 */
void
TC_REP::do_matrix_assignment (const tree_constant& rhs,
			      Range& ri, idx_vector& jv)
{
  REP_RHS_MATRIX (rhs, rhs_m, rhs_cm, rhs_nr, rhs_nc);

  double b = ri.base ();
  double increment = ri.inc ();

  for (int j = 0; j < jv.capacity (); j++)
    {
      int col = jv.elem (j);
      for (int i = 0; i < ri.nelem (); i++)
	{
	  double tmp = b + i * increment;
	  int row = tree_to_mat_idx (tmp);
	  REP_ELEM_ASSIGN (row, col, rhs_m.elem (i, j),
			   rhs_m.elem (i, j), rhs.is_real_type ());
	}
    }
}

/* 11 */
void
TC_REP::do_matrix_assignment (const tree_constant& rhs,
			      Range& ri, Range& rj)
{
  double ib = ri.base ();
  double iinc = ri.inc ();
  double jb = rj.base ();
  double jinc = rj.inc ();

  REP_RHS_MATRIX (rhs, rhs_m, rhs_cm, rhs_nr, rhs_nc);

  for (int i = 0; i < ri.nelem (); i++)
    {
      double itmp = ib + i * iinc;
      int row = tree_to_mat_idx (itmp);
      for (int j = 0; j < rj.nelem (); j++)
	{
	  double jtmp = jb + j * jinc;
	  int col = tree_to_mat_idx (jtmp);
	  REP_ELEM_ASSIGN (row, col, rhs_m.elem  (i, j),
			   rhs_cm.elem (i, j), rhs.is_real_type ());
	}
    }
}

/* 12 */
void
TC_REP::do_matrix_assignment (const tree_constant& rhs,
			      Range& ri, TC_REP::constant_type mcj)
{
  assert (mcj == magic_colon);

  if (rhs.is_zero_by_zero ())
    {
      delete_rows (ri);
    }
  else
    {
      REP_RHS_MATRIX (rhs, rhs_m, rhs_cm, rhs_nr, rhs_nc);

      double ib = ri.base ();
      double iinc = ri.inc ();

      int nc = columns ();

      for (int i = 0; i < ri.nelem (); i++)
	{
	  double itmp = ib + i * iinc;
	  int row = tree_to_mat_idx (itmp);
	  for (int j = 0; j < nc; j++)
	    REP_ELEM_ASSIGN (row, j, rhs_m.elem (i, j),
			     rhs_cm.elem (i, j), rhs.is_real_type ());
	}
    }
}

/* 13 */
void
TC_REP::do_matrix_assignment (const tree_constant& rhs,
			      TC_REP::constant_type mci, int j)
{
  assert (mci == magic_colon);

  int nr = rows ();

  if (rhs.is_zero_by_zero ())
    {
      delete_column (j);
    }
  else if (rhs.is_matrix_type ())
    {
      REP_RHS_MATRIX (rhs, rhs_m, rhs_cm, rhs_nr, rhs_nc);

      for (int i = 0; i < nr; i++)
	REP_ELEM_ASSIGN (i, j, rhs_m.elem (i, 0),
			 rhs_cm.elem (i, 0), rhs.is_real_type ());
    }
  else if (rhs.is_scalar_type () && nr == 1)
    {
      REP_ELEM_ASSIGN (0, j, rhs.double_value (),
		       rhs.complex_value (), rhs.is_real_type ());
    }
  else
    panic_impossible ();
}

/* 14 */
void
TC_REP::do_matrix_assignment (const tree_constant& rhs,
			      TC_REP::constant_type mci, idx_vector& jv)
{
  assert (mci == magic_colon);

  if (rhs.is_zero_by_zero ())
    {
      delete_columns (jv);
    }
  else
    {
      REP_RHS_MATRIX (rhs, rhs_m, rhs_cm, rhs_nr, rhs_nc);

      int nr = rows ();

      for (int i = 0; i < nr; i++)
	{
	  for (int j = 0; j < jv.capacity (); j++)
	    {
	      int col = jv.elem (j);
	      REP_ELEM_ASSIGN (i, col, rhs_m.elem (i, j),
			       rhs_cm.elem (i, j), rhs.is_real_type ());
	    }
	}
    }
}

/* 15 */
void
TC_REP::do_matrix_assignment (const tree_constant& rhs,
			      TC_REP::constant_type mci, Range& rj)
{
  assert (mci == magic_colon);

  if (rhs.is_zero_by_zero ())
    {
      delete_columns (rj);
    }
  else
    {
      REP_RHS_MATRIX (rhs, rhs_m, rhs_cm, rhs_nr, rhs_nc);

      int nr = rows ();

      double jb = rj.base ();
      double jinc = rj.inc ();

      for (int j = 0; j < rj.nelem (); j++)
	{
	  double jtmp = jb + j * jinc;
	  int col = tree_to_mat_idx (jtmp);
	  for (int i = 0; i < nr; i++)
	    {
	      REP_ELEM_ASSIGN (i, col, rhs_m.elem (i, j),
			       rhs_cm.elem (i, j), rhs.is_real_type ());
	    }
	}
    }
}

/* 16 */
void
TC_REP::do_matrix_assignment (const tree_constant& rhs,
			      TC_REP::constant_type mci,
			      TC_REP::constant_type mcj)
{
  assert (mci == magic_colon && mcj == magic_colon);

  switch (type_tag)
    {
    case scalar_constant:
      break;

    case matrix_constant:
      delete matrix;
      break;

    case complex_scalar_constant:
      delete complex_scalar;
      break;

    case complex_matrix_constant:
      delete complex_matrix;
      break;

    case string_constant:
      delete [] string;
      break;

    case range_constant:
      delete range;
      break;

    case magic_colon:
    default:
      panic_impossible ();
      break;
    }

  type_tag = rhs.const_type ();

  switch (type_tag)
    {
    case scalar_constant:
      scalar = rhs.double_value ();
      break;

    case matrix_constant:
      matrix = new Matrix (rhs.matrix_value ());
      break;

    case string_constant:
      string = strsave (rhs.string_value ());
      break;

    case complex_matrix_constant:
      complex_matrix = new ComplexMatrix (rhs.complex_matrix_value ());
      break;

    case complex_scalar_constant:
      complex_scalar = new Complex (rhs.complex_value ());
      break;

    case range_constant:
      range = new Range (rhs.range_value ());
      break;

    case magic_colon:
    default:
      panic_impossible ();
      break;
    }
}

// Functions for deleting rows or columns of a matrix.  These are used
// to handle statements like
//
//   M (i, j) = []

void
TC_REP::delete_row (int idx)
{
  if (type_tag == matrix_constant)
    {
      int nr = matrix->rows ();
      int nc = matrix->columns ();
      Matrix *new_matrix = new Matrix (nr-1, nc);
      int ii = 0;
      for (int i = 0; i < nr; i++)
	{
	  if (i != idx)
	    {
	      for (int j = 0; j < nc; j++)
		new_matrix->elem (ii, j) = matrix->elem (i, j);
	      ii++;
	    }
	}
      delete matrix;
      matrix = new_matrix;
    }
  else if (type_tag == complex_matrix_constant)
    {
      int nr = complex_matrix->rows ();
      int nc = complex_matrix->columns ();
      ComplexMatrix *new_matrix = new ComplexMatrix (nr-1, nc);
      int ii = 0;
      for (int i = 0; i < nr; i++)
	{
	  if (i != idx)
	    {
	      for (int j = 0; j < nc; j++)
		new_matrix->elem (ii, j) = complex_matrix->elem (i, j);
	      ii++;
	    }
	}
      delete complex_matrix;
      complex_matrix = new_matrix;
    }
  else
    panic_impossible ();
}

void
TC_REP::delete_rows (idx_vector& iv)
{
  iv.sort_uniq ();
  int num_to_delete = iv.length ();

  if (num_to_delete == 0)
    return;

  int nr = rows ();
  int nc = columns ();

// If deleting all rows of a column vector, make result 0x0.
  if (nc == 1 && num_to_delete == nr)
    nc = 0;

  if (type_tag == matrix_constant)
    {
      Matrix *new_matrix = new Matrix (nr-num_to_delete, nc);
      if (nr > num_to_delete)
	{
	  int ii = 0;
	  int idx = 0;
	  for (int i = 0; i < nr; i++)
	    {
	      if (i == iv.elem (idx))
		idx++;
	      else
		{
		  for (int j = 0; j < nc; j++)
		    new_matrix->elem (ii, j) = matrix->elem (i, j);
		  ii++;
		}
	    }
	}
      delete matrix;
      matrix = new_matrix;
    }
  else if (type_tag == complex_matrix_constant)
    {
      ComplexMatrix *new_matrix = new ComplexMatrix (nr-num_to_delete, nc);
      if (nr > num_to_delete)
	{
	  int ii = 0;
	  int idx = 0;
	  for (int i = 0; i < nr; i++)
	    {
	      if (i == iv.elem (idx))
		idx++;
	      else
		{
		  for (int j = 0; j < nc; j++)
		    new_matrix->elem (ii, j) = complex_matrix->elem (i, j);
		  ii++;
		}
	    }
	}
      delete complex_matrix;
      complex_matrix = new_matrix;
    }
  else
    panic_impossible ();
}

void
TC_REP::delete_rows (Range& ri)
{
  ri.sort ();
  int num_to_delete = ri.nelem ();

  if (num_to_delete == 0)
    return;

  int nr = rows ();
  int nc = columns ();

// If deleting all rows of a column vector, make result 0x0.
  if (nc == 1 && num_to_delete == nr)
    nc = 0;

  double ib = ri.base ();
  double iinc = ri.inc ();

  int max_idx = tree_to_mat_idx (ri.max ());

  if (type_tag == matrix_constant)
    {
      Matrix *new_matrix = new Matrix (nr-num_to_delete, nc);
      if (nr > num_to_delete)
	{
	  int ii = 0;
	  int idx = 0;
	  for (int i = 0; i < nr; i++)
	    {
	      double itmp = ib + idx * iinc;
	      int row = tree_to_mat_idx (itmp);

	      if (i == row && row <= max_idx)
		idx++;
	      else
		{
		  for (int j = 0; j < nc; j++)
		    new_matrix->elem (ii, j) = matrix->elem (i, j);
		  ii++;
		}
	    }
	}
      delete matrix;
      matrix = new_matrix;
    }
  else if (type_tag == complex_matrix_constant)
    {
      ComplexMatrix *new_matrix = new ComplexMatrix (nr-num_to_delete, nc);
      if (nr > num_to_delete)
	{
	  int ii = 0;
	  int idx = 0;
	  for (int i = 0; i < nr; i++)
	    {
	      double itmp = ib + idx * iinc;
	      int row = tree_to_mat_idx (itmp);

	      if (i == row && row <= max_idx)
		idx++;
	      else
		{
		  for (int j = 0; j < nc; j++)
		    new_matrix->elem (ii, j) = complex_matrix->elem (i, j);
		  ii++;
		}
	    }
	}
      delete complex_matrix;
      complex_matrix = new_matrix;
    }
  else
    panic_impossible ();
}

void
TC_REP::delete_column (int idx)
{
  if (type_tag == matrix_constant)
    {
      int nr = matrix->rows ();
      int nc = matrix->columns ();
      Matrix *new_matrix = new Matrix (nr, nc-1);
      int jj = 0;
      for (int j = 0; j < nc; j++)
	{
	  if (j != idx)
	    {
	      for (int i = 0; i < nr; i++)
		new_matrix->elem (i, jj) = matrix->elem (i, j);
	      jj++;
	    }
	}
      delete matrix;
      matrix = new_matrix;
    }
  else if (type_tag == complex_matrix_constant)
    {
      int nr = complex_matrix->rows ();
      int nc = complex_matrix->columns ();
      ComplexMatrix *new_matrix = new ComplexMatrix (nr, nc-1);
      int jj = 0;
      for (int j = 0; j < nc; j++)
	{
	  if (j != idx)
	    {
	      for (int i = 0; i < nr; i++)
		new_matrix->elem (i, jj) = complex_matrix->elem (i, j);
	      jj++;
	    }
	}
      delete complex_matrix;
      complex_matrix = new_matrix;
    }
  else
    panic_impossible ();
}

void
TC_REP::delete_columns (idx_vector& jv)
{
  jv.sort_uniq ();
  int num_to_delete = jv.length ();

  if (num_to_delete == 0)
    return;

  int nr = rows ();
  int nc = columns ();

// If deleting all columns of a row vector, make result 0x0.
  if (nr == 1 && num_to_delete == nc)
    nr = 0;

  if (type_tag == matrix_constant)
    {
      Matrix *new_matrix = new Matrix (nr, nc-num_to_delete);
      if (nc > num_to_delete)
	{
	  int jj = 0;
	  int idx = 0;
	  for (int j = 0; j < nc; j++)
	    {
	      if (j == jv.elem (idx))
		idx++;
	      else
		{
		  for (int i = 0; i < nr; i++)
		    new_matrix->elem (i, jj) = matrix->elem (i, j);
		  jj++;
		}
	    }
	}
      delete matrix;
      matrix = new_matrix;
    }
  else if (type_tag == complex_matrix_constant)
    {
      ComplexMatrix *new_matrix = new ComplexMatrix (nr, nc-num_to_delete);
      if (nc > num_to_delete)
	{
	  int jj = 0;
	  int idx = 0;
	  for (int j = 0; j < nc; j++)
	    {
	      if (j == jv.elem (idx))
		idx++;
	      else
		{
		  for (int i = 0; i < nr; i++)
		    new_matrix->elem (i, jj) = complex_matrix->elem (i, j);
		  jj++;
		}
	    }
	}
      delete complex_matrix;
      complex_matrix = new_matrix;
    }
  else
    panic_impossible ();
}

void
TC_REP::delete_columns (Range& rj)
{
  rj.sort ();
  int num_to_delete = rj.nelem ();

  if (num_to_delete == 0)
    return;

  int nr = rows ();
  int nc = columns ();

// If deleting all columns of a row vector, make result 0x0.
  if (nr == 1 && num_to_delete == nc)
    nr = 0;

  double jb = rj.base ();
  double jinc = rj.inc ();

  int max_idx = tree_to_mat_idx (rj.max ());

  if (type_tag == matrix_constant)
    {
      Matrix *new_matrix = new Matrix (nr, nc-num_to_delete);
      if (nc > num_to_delete)
	{
	  int jj = 0;
	  int idx = 0;
	  for (int j = 0; j < nc; j++)
	    {
	      double jtmp = jb + idx * jinc;
	      int col = tree_to_mat_idx (jtmp);

	      if (j == col && col <= max_idx)
		idx++;
	      else
		{
		  for (int i = 0; i < nr; i++)
		    new_matrix->elem (i, jj) = matrix->elem (i, j);
		  jj++;
		}
	    }
	}
      delete matrix;
      matrix = new_matrix;
    }
  else if (type_tag == complex_matrix_constant)
    {
      ComplexMatrix *new_matrix = new ComplexMatrix (nr, nc-num_to_delete);
      if (nc > num_to_delete)
	{
	  int jj = 0;
	  int idx = 0;
	  for (int j = 0; j < nc; j++)
	    {
	      double jtmp = jb + idx * jinc;
	      int col = tree_to_mat_idx (jtmp);

	      if (j == col && col <= max_idx)
		idx++;
	      else
		{
		  for (int i = 0; i < nr; i++)
		    new_matrix->elem (i, jj) = complex_matrix->elem (i, j);
		  jj++;
		}
	    }
	}
      delete complex_matrix;
      complex_matrix = new_matrix;
    }
  else
    panic_impossible ();
}

/*
;;; Local Variables: ***
;;; mode: C++ ***
;;; page-delimiter: "^/\\*" ***
;;; End: ***
*/