Index: mttroot/mtt/bin/mtt
==================================================================
--- mttroot/mtt/bin/mtt
+++ mttroot/mtt/bin/mtt
@@ -13,10 +13,13 @@
 ###############################################################
 ## Version control history
 ###############################################################
 ## $Header$
 ## $Log$
+## Revision 1.292  2001/02/05 17:27:40  gawthrop
+## Make sure _def.r exists before creating _state.txt
+##
 ## Revision 1.291  2000/12/27 14:50:40  peterg
 ## This is the first CVS version (4.9).
 ## Commented out code now deleted
 ##
 ## Revision 1.290  2000/12/05 09:59:37  peterg
@@ -1747,13 +1750,13 @@
     
    if [ -f "$filename" ]; then
       echo $filename exists
    else
      if [ -n "$Verbose" ]; then
-      echo make -s -f $2_rep.make "SYS=$1" "LANG=$3" "ARG=$4" 
+      echo make -s -f $2_rep.make "SYS=$1" "LANG=$3" "ARG=$4" "OPTS=$mtt_switches" 
      fi
-      make -s -f $2_rep.make "SYS=$1" "LANG=$3" "ARG=$4"
+      make -s -f $2_rep.make "SYS=$1" "LANG=$3" "ARG=$4" "OPTS=$mtt_switches"
       if [ -n "$4" ]; then
          echo Copying  $1_$2$__ARGS.$ps
          cp  $1_$2$__ARGS.$ps ..
       fi
    fi

ADDED   mttroot/mtt/lib/cc/mtt_euler.cc
Index: mttroot/mtt/lib/cc/mtt_euler.cc
==================================================================
--- /dev/null
+++ mttroot/mtt/lib/cc/mtt_euler.cc
@@ -0,0 +1,36 @@
+#include <octave/oct.h>
+
+DEFUN_DLD (mtt_euler, args, ,
+	   "euler integration method")
+{
+#ifdef OCTAVE_DEV
+  ColumnVector		x	= args(0).column_vector_value ();
+  const ColumnVector	dx	= args(1).column_vector_value ();
+  const double		ddt	= args(2).double_value ();
+  const int		Nx	= static_cast<int> (args(3).double_value ());
+  const ColumnVector   	openx	= args(4).column_vector_value ();
+#else
+  const ColumnVector	x	= args(0).vector_value ();
+  const ColumnVector   	dx	= args(1).vector_value ();
+  const double		ddt	= args(2).double_value ();
+  const int		Nx	= static_cast<int> (args(3).double_value ());
+  const ColumnVector   	openx	= args(4).vector_value ();
+#endif
+
+  register int i, n;
+  
+  n = Nx;
+  for (i = 0; i < Nx; i++)
+    {
+      if (0 != openx (i))
+	{
+	  x (i) = 0.0;
+	}
+      else
+	{
+	  x (i) += dx (i) * ddt;
+	}
+    }
+  
+  return octave_value (x);
+}

ADDED   mttroot/mtt/lib/cc/mtt_implicit.cc
Index: mttroot/mtt/lib/cc/mtt_implicit.cc
==================================================================
--- /dev/null
+++ mttroot/mtt/lib/cc/mtt_implicit.cc
@@ -0,0 +1,139 @@
+#include <octave/oct.h>
+#include "useful-functions.hh"
+
+#include <octave/xdiv.h>
+static inline int
+result_ok (int info, double rcond, int warn = 1)
+{
+  assert (info != -1);
+
+  if (info == -2)
+    {
+      if (warn)
+	warning ("matrix singular to machine precision, rcond = %g", rcond);
+      else
+	error ("matrix singular to machine precision, rcond = %g", rcond);
+
+      return 0;
+    }
+  else
+    return 1;
+}
+bool
+mx_leftdiv_conform (const Matrix& a, const ColumnVector& b)
+{
+  int a_nr = a.rows ();
+  int b_nr = b.length ();
+
+  if (a_nr != b_nr)
+    {
+      int a_nc = a.cols ();
+      int b_nc = 1;
+
+      gripe_nonconformant ("operator \\", a_nr, a_nc, b_nr, b_nc);
+      return false;
+    }
+
+  return true;
+}
+
+// need to update xdiv.cc ?
+inline ColumnVector
+ldiv (const Matrix &a, const ColumnVector &b)
+{
+  if (! mx_leftdiv_conform (a, b))
+    return ColumnVector ();
+
+  int info;
+  if (a.rows () == a.columns ())
+    {
+      double rcond = 0.0;
+      ColumnVector result = a.solve (b, info, rcond);
+      if (result_ok (info, rcond))
+	return result;
+    }
+
+  int rank;
+  return a.lssolve (b, info, rank);
+}
+
+DEFUN_DLD (mtt_implicit, args, ,
+	   "implicit integration method")
+{
+#ifdef OCTAVE_DEV
+  ColumnVector		x	= args(0).column_vector_value ();
+  const ColumnVector	dx	= args(1).column_vector_value ();
+  const Matrix		AA	= args(2).matrix_value ();
+  const ColumnVector	AAx	= args(3).column_vector_value ();
+  const double		t	= args(4).double_value ();
+  const int		Nx	= (int) (args(5).double_value ());
+  const ColumnVector	openx	= args(6).column_vector_value ();
+#else
+  ColumnVector		x	= args(0).vector_value ();
+  const ColumnVector	dx	= args(1).vector_value ();
+  const Matrix		AA	= args(2).matrix_value ();
+  const ColumnVector	AAx	= args(3).vector_value ();
+  const double		t	= args(4).double_value ();
+  const int		Nx	= (int) (args(5).double_value ());
+  const ColumnVector	openx	= args(6).vector_value ();
+#endif
+
+  register int i, n;
+  register int col_old, col_new;
+  register int row_old, row_new;
+
+  n = Nx;
+  for (i = 0; i < Nx; i++)
+    {
+      if (0 != openx (i))
+	{
+	  n--;
+	}
+    }
+
+  ColumnVector	tmp_dx	(n, 0.0);
+  ColumnVector	tmp_x	(n, 0.0);
+  ColumnVector	tmp_AAx	(n, 0.0);
+  Matrix	tmp_AA	(n, n, 0.0);
+
+  for (row_new = row_old = 0; row_old < Nx; row_old++)
+    {
+      if (0 == openx (row_old))
+	{
+	  tmp_dx  (row_new)	= dx  (row_old);
+	  tmp_AAx (row_new)	= AAx (row_old);
+	  for (col_new = col_old = 0; col_old < Nx; col_old++)
+	    {
+	      if (0 == openx (col_old))
+		{
+		  // xxx: this can be improved by symmetry
+		  tmp_AA (row_new,col_new) = AA (row_old,col_old);
+		  col_new++;
+		}
+	    }
+	  row_new++;
+	}
+    }
+
+  // can't get ldiv to work - doesn't like ColVector
+  // tmp_x = tmp_AA.pseudo_inverse () * tmp_AAx + t * tmp_dx;
+  tmp_x = ldiv (tmp_AA, (tmp_AAx + t * tmp_dx));
+  
+  row_new = 0;
+  for (row_old = 0; row_old < Nx; row_old++)
+    {
+      if (0 == openx (row_old))
+	{
+	  x (row_old) = tmp_x (row_new);
+	  row_new++;
+	}
+      else
+	{
+	  x (row_old) = 0.0;
+	}
+    }
+
+  return octave_value (x);
+}
+
+

ADDED   mttroot/mtt/lib/cc/mtt_write.cc
Index: mttroot/mtt/lib/cc/mtt_write.cc
==================================================================
--- /dev/null
+++ mttroot/mtt/lib/cc/mtt_write.cc
@@ -0,0 +1,54 @@
+#include <octave/oct.h>
+#include <octave/variables.h>
+
+#ifdef STANDALONE
+static Matrix MTT_data;
+#endif
+
+DEFUN_DLD (mtt_write, args, ,
+	   "append current data to output")
+{
+  const double		t	= args(0).double_value ();
+#ifdef OCTAVE_DEV
+  const ColumnVector	x	= args(1).column_vector_value ();
+  const ColumnVector	y	= args(2).column_vector_value ();
+#else
+  const ColumnVector	x	= args(1).vvector_value ();
+  const ColumnVector	y	= args(2).vector_value ();
+#endif
+
+  const int		nx	= x.length ();
+  const int		ny	= y.length ();
+
+
+  ColumnVector Output (2+nx+ny, 0.0);
+  
+  Output (0) = Output (1+nx) = t;
+  Output.insert (x.transpose (), 1);
+  Output.insert (y.transpose (), 2+nx);
+
+  Matrix data;
+
+  if (0.0 == t)
+    {
+      data = static_cast<Matrix> (Output.transpose ());
+    }
+  else
+    {
+#ifdef STANDALONE
+      data = MTT_data.transpose ();
+#else
+      data = get_global_value ("MTT_data").matrix_value ().transpose ();
+#endif
+      data = data.append (Output);
+    }
+  data = data.transpose ();
+#ifdef STANDALONE
+  MTT_data = data;
+  cout << Output.transpose () << endl;
+#else
+  set_global_value ("MTT_data", data);
+#endif
+  return data;
+}
+

ADDED   mttroot/mtt/lib/rep/standalone_rep.make
Index: mttroot/mtt/lib/rep/standalone_rep.make
==================================================================
--- /dev/null
+++ mttroot/mtt/lib/rep/standalone_rep.make
@@ -0,0 +1,57 @@
+# -*-makefile-*-
+
+.POSIX:
+
+MTTFLAGS	= -q -u -oct $(OPTS)
+
+# Adapt according to local set-up and mkoctfile
+CC		= g++
+CXXFLAGS	= $(DEBUG) $(OPTIM) $(DEFINES) $(ARCHFLAGS) -fno-rtti -fno-exceptions -fno-implicit-templates
+
+DEBUG		= -v -pg -g
+OPTIM		= -O3
+
+PREFIX		= /usr/local
+
+INCLUDES	= -I$(PREFIX)/include/octave
+
+OCTAVE_SRC_PATH	= $(PREFIX)/src/octave
+
+LIBOCTAVE	= -L$(PREFIX)/lib/octave			-loctave -lcruft -loctinterp
+LIBKPATHSEA	= -L$(OCTAVE_SRC_PATH)/kpathsea			-lkpathsea
+LIBREADLINE	= -L$(OCTAVE_SRC_PATH)/readline			-lreadline
+LIBBLAS		= -L/usr/local/src/ATLAS/lib/Linux_PIII		-lcblas -lf77blas -llapack -latlas -ltstatlas
+LIBF2C		= 						-lg2c
+LIBRARIES	= 						-ldl -lm -lncurses
+
+ARCHFLAGS	= $(i386FLAGS)
+i386FLAGS	=  -mieee-fp 
+
+# Define -DOCTAVE_DEV for octave 2.1.x
+ifeq (0, $(shell octave --version | awk -F\. '{print $2}'))
+DEFINES   = -DSTANDALONE
+else
+DEFINES   = -DSTANDALONE -DOCTAVE_DEV
+endif
+
+all: $(SYS)_standalone.$(LANG)
+
+$(SYS)_standalone.exe: $(SYS)_ode2odes.cc $(SYS)_def.h $(SYS)_sympar.h
+	cp $(MTT_LIB)/cc/*.cc .
+	$(CC) *.cc -o $@ $(CXXFLAGS) $(INCLUDES) $(LIBOCTAVE) $(LIBKPATHSEA) $(LIBREADLINE) $(LIBBLAS) $(LIBF2C) $(LIBRARIES)
+
+.PHONY: $(SYS)_standalone.clean
+
+$(SYS)_ode2odes.cc:
+	mtt $(MTTFLAGS) $(SYS) ode2odes m
+
+$(SYS)_def.h:
+	mtt $(MTTFLAGS) $(SYS) def h
+
+$(SYS)_sympar.h:
+	mtt $(MTTFLAGS) $(SYS) sympar h 
+
+
+$(SYS)_standalone.clean:
+	cd .. ; mtt Clean
+	rm -f ../$(SYS)_standalone.exe