# Copyright (C) 2000 by Peter J. Gawthrop

###############################################################
## Version control history
###############################################################
## $Header$
## $Log$
## Revision 1.293.2.2  2001/02/23 03:53:53  geraint
## Minor revisions. (ode2odes.exe works on GNU/Linux and ~GNU/Irix)
##
## Revision 1.293.2.1  2001/02/22 06:41:40  geraint
## Initial merge of ode2odes.exe into main mtt.
## standalone_rep.make deleted: rules moved to mtt, variables to mttrc.
##
## Revision 1.293  2001/02/06 04:39:08  geraint
## Added _standalone.exe representation.
##

%.dvi: %.tex

# MTT implicit rules
## .oct files
$1_%.oct:  $1_%.cc $1_def.h $1_sympar.h $1_cr.h
	echo Creating $1_\$*.oct; $MKOCTFILE  $1_\$*.cc

.PRECIOUS: mtt_%.oct
mtt_%.oct:: mtt_%.cc
	echo Compiling \$<
	${MKOCTFILE} $define_octave_dev \$<

## .cc files
.PRECIOUS: %.cc # Don't let mtt delete them
$1_%.cc:  $1_%.m
	mtt_m2cc.sh  $1 \$* cat 

mtt_%.cc: ${MTT_LIB}/cc/mtt_%.cc
	cp ${MTT_LIB}/cc/mtt_\$*.cc mtt_\$*.cc

## .o files
.PRECIOUS: $1_%.o
$1_%.o: $1_%.cc $1_def.h $1_sympar.h $1_cr.h
	echo Compiling $1_\$*.cc
	${MTT_CXX} ${MTT_CXXFLAGS} ${MTT_CXXINCS} -c $< -DSTANDALONE

mtt_%.o:: ${MTT_LIB}/cc/mtt_%.cc
mtt_%.o: mtt_%.cc
	echo Compiling mtt_\$*.cc
	${MTT_CXX} ${MTT_CXXFLAGS} ${MTT_CXXINCS} -c $< -DSTANDALONE $define_octave_dev

# Generates code to write reduce code
$1_ode_write.r: $1_def.m
	def2write_r  $1 $2

#SUMMARY simpar	simulation information (m)
$1_simpar.m: $1_simpar.txt $1_def.r $1_sympar.txt
	mtt_txt2m $1 simpar


#SUMMARY numpar	numerical parameter declaration -- default (txt) 
$1_numpar.txt : FORCE
$1_numpar.txt : $1_sympar.txt $1_def.r
	mtt -q $1 sympar txt
	mtt -q $1 def r
ifeq ($target,$1_numpar.txt)
	mtt_update $1 numpar update
else
	mtt_update $1 numpar
endif

#SUMMARY state	numerical state declaration  -- default (txt)
$1_state.txt : FORCE
$1_state.txt : $1_struc.txt $1_sympar.txt $1_def.r
	mtt -q $1 struc txt
	mtt -q $1 sympar txt
	mtt -q $1 def r
ifeq ($target,$1_state.txt)
	mtt_update $1 state update
else
	mtt_update $1 state
endif

#SUMMARY input	numerical state declaration  -- default (txt)
$1_input.txt : FORCE
$1_input.txt : $1_struc.txt $1_sympar.txt
	mtt -q $1 struc txt
	mtt -q $1 sympar txt
ifeq ($target,$1_input.txt)
	mtt_update $1 input update
else
	mtt_update $1 input
endif

#SUMMARY logic	Dynamic switch logic (txt)
$1_logic.txt : FORCE
$1_logic.txt : $1_def.r $1_struc.txt $1_sympar.txt
	mtt -q $1 def r
	mtt -q $1 struc txt
	mtt -q $1 sympar txt
ifeq ($target,$1_logic.txt)
	mtt_update $1 logic update
else
	mtt_update $1 logic
endif

# Dummy target
FORCE:

ode2odes_common.m    : $1_input.m   $1_logic.m   $1_numpar.m   $1_simpar.m   $1_state.m                              ; @echo > /dev/null
ode2odes_common.cc   : $1_input.cc  $1_logic.cc  $1_numpar.cc  $1_simpar.cc  $1_state.cc                             ; @echo > /dev/null
ode2odes_common.o    : $1_input.o   $1_logic.o   $1_numpar.o   $1_simpar.o   $1_state.o   mtt_euler.o mtt_implicit.o ; @echo > /dev/null
ode2odes_common.oct  : $1_input.oct $1_logic.oct $1_numpar.oct $1_simpar.oct $1_state.oct                            ; @echo > /dev/null
ode2odes_euler.m     : $1_ode.m     $1_odeo.m                                                                        ; @echo > /dev/null
ode2odes_euler.cc    : $1_ode.cc    $1_odeo.cc                                                                       ; @echo > /dev/null
ode2odes_euler.o     : $1_ode.o     $1_odeo.o                                                                        ; @echo > /dev/null
ode2odes_euler.oct   : $1_ode.oct   $1_odeo.oct                                                                      ; @echo > /dev/null
ode2odes_implicit.m  : $1_cseo.m    $1_csex.m   $1_smxa.m   $1_smxax.m                                               ; @echo > /dev/null
ode2odes_implicit.cc : $1_cseo.cc   $1_csex.cc  $1_smxa.cc  $1_smxax.cc                                              ; @echo > /dev/null
ode2odes_implicit.o  : $1_cseo.o    $1_csex.o   $1_smxa.o   $1_smxax.o                                               ; @echo > /dev/null
ode2odes_implicit.oct: $1_cseo.oct  $1_csex.oct $1_smxa.oct $1_smxax.oct                                             ; @echo > /dev/null
ode2odes_common.m    : $1_input.m   $1_logic.m   $1_numpar.m   $1_simpar.m   $1_state.m
	@echo > /dev/null
ode2odes_common.cc   : $1_input.cc  $1_logic.cc  $1_numpar.cc  $1_simpar.cc  $1_state.cc
	@echo > /dev/null
ode2odes_common.o    : $1_input.o   $1_logic.o   $1_numpar.o   $1_simpar.o   $1_state.o   mtt_euler.o   mtt_implicit.o
	@echo > /dev/null
ode2odes_common.oct  : $1_input.oct $1_logic.oct $1_numpar.oct $1_simpar.oct $1_state.oct mtt_euler.oct mtt_implicit.oct
	@echo > /dev/null
ode2odes_euler.m     : $1_ode.m     $1_odeo.m
	@echo > /dev/null
ode2odes_euler.cc    : $1_ode.cc    $1_odeo.cc
	@echo > /dev/null
ode2odes_euler.o     : $1_ode.o     $1_odeo.o
	@echo > /dev/null
ode2odes_euler.oct   : $1_ode.oct   $1_odeo.oct
	@echo > /dev/null
ode2odes_implicit.m  : $1_cseo.m    $1_csex.m   $1_smxa.m   $1_smxax.m
	@echo > /dev/null
ode2odes_implicit.cc : $1_cseo.cc   $1_csex.cc  $1_smxa.cc  $1_smxax.cc
	@echo > /dev/null
ode2odes_implicit.o  : $1_cseo.o    $1_csex.o   $1_smxa.o   $1_smxax.o
	@echo > /dev/null
ode2odes_implicit.oct: $1_cseo.oct  $1_csex.oct $1_smxa.oct $1_smxax.oct
	@echo > /dev/null


#SUMMARY numpar	numerical parameter declaration (m) 
$1_numpar.m:  $1_numpar.txt $1_sympars.txt
	mtt_txt2m $1 numpar

#SUMMARY numpar	numerical parameter declaration (c) 

endif
endif
ifeq ($integration_method,euler)
$1_ode2odes.m : $1_def.r $1_sympars.txt\
		$1_simpar.m $1_numpar.m $1_state.m $1_input.m \
		$1_ode.m $1_odeo.m  $1_logic.m
ifeq ($using_oct,yes)
	echo "*** Warning: Shouldn't be here! Creating dummy $1_ode2odes.m"
	touch $1_ode2odes.m # Create a dummy which wont' be used
	mtt $mtt_switches -q -u $1 ode2odes oct
else
	make_ode2odes $1 m $integration_method
endif

endif

	m2dat $1_sro
1_odes.c: $1_ode.c $1_def.r $1_sympar.r
	ode2odes_r2c $1 
$1_sms.c: $1_sm.c $1_def.r $1_sympar.r
	sm2sms_r2c $1

ifeq ($computation,octave)
ifeq ($using_oct,yes)
$1_odes.dat2: $1_ode2odes.oct
	time octave_ode2odes $1
else
$1_odes.dat2: $1_ode2odes.m
	time octave_ode2odes $1
endif
endif

ifeq ($computation,c)
$1_ode2odes.out: $1_ode2odes.c
	echo Creating $1_odes.out
	$CC $1_ode2odes.c -lp2c -lm -o $1_ode2odes.out

###############################################################
## Version control history
###############################################################
## $Id$
## $Log$
## Revision 1.17.2.2  2001/02/23 03:53:53  geraint
## Minor revisions. (ode2odes.exe works on GNU/Linux and ~GNU/Irix)
##
## Revision 1.17.2.1  2001/02/22 06:41:40  geraint
## Initial merge of ode2odes.exe into main mtt.
## standalone_rep.make deleted: rules moved to mtt, variables to mttrc.
##
## Revision 1.17  2000/12/27 16:46:13  peterg
## Stripped the mtt- from paths
##

    # exported variables

    export MTT_CXX="g++"
    export MTT_CXXFLAGS="${DEBUG} ${OPTIM} ${FLAGS}"
    export MTT_CXXLIBS="${LOCTAVE} ${LKPATHSEA} ${LREADLINE} ${LF2C} ${LSYSTEM}"
    export MTT_CXXINCS="-I. ${IOCTAVE}"

    export MTT_LDFLAGS=""
fi

#! /bin/sh

     ###################################### 
     ##### Model Transformation Tools #####
    ######################################

###############################################################
## Version control history
###############################################################
## $Id$
## $Log$
## Revision 1.51  2001/02/19 06:33:19  geraint
## Removed operation form loop.
##
## Revision 1.50  2001/02/18 09:18:49  geraint
## Removed temporary Matrices from mtt_implicit.cc
##
## Revision 1.49  2001/02/14 06:06:34  geraint
## Removed octave_value_list wrappers from standalone.exe - speed improvements
##
## Revision 1.48  2001/02/11 07:08:59  geraint

    vector_value=vector_value
    feval_header=toplev.h
fi

cat <<EOF  > $filename
#include <octave/oct.h>
#include <octave/ov-struct.h>
#include <octave/oct-map.h>
#include <octave/load-save.h>
#include <octave/lo-mappers.h>
#include <octave/variables.h>

#ifndef STANDALONE
#include <octave/${feval_header}>
#endif

  data.elem (row, ++col) = t;
  for (register int i = 0; i < nx; i++)
    data.elem (row, ++col) = x.elem (i);

  row++;

  if (nrows == row)
    {
#ifdef STANDALONE
      save_ascii_data_for_plotting (cout, data, "MTT_data");
    cout << data << endl;
//    cout << data << endl;
#else // ! STANDALONE
  set_global_value ("MTT_data", data);
#endif

    }
}

#ifdef STANDALONE
int main (void) {
#else
DEFUN_DLD (${sys}_ode2odes, args, ,
"Octave ode2odes representation of system 

  static ColumnVector	open_switches (MTTNX);

  register double t	= 0.0;

  const double	ddt	= dt / stepfactor;
  const int	ilast	= static_cast<int> (round ((last - first) / ddt)) + 1;
  const int	nrows	= static_cast<int> (round ((last - first) / dt)) + 1;

  for (register int i = 0; i < MTTNU; i++)
    {
      u (i) = 0.0;
    }

  for (register int j = 0, i = 1; i <= ilast; i++)
    {
      y	= mtt_${odeo} (x, u, t, par);
      u	= mtt_input (x, y, t, par);
      if (0 == j)
	{

#ifndef HAVE_USEFUL_FUNCTIONS_HH
#define HAVE_USEFUL_FUNCTIONS_HH

template <class T>
inline T max (const T &x1, const T &x2)
{
  return (x1 >= x2) ? x1 : (x2 < x1) ? x2 : 0;
}

sign (T x)
{
  return
    (0 < x) ? +1 :
    (0 > x) ? -1 :
    0;
}

#endif // HAVE_USEFUL_FUNCTIONS_HH

{
#else // !STANDALONE
DEFUN_DLD (mtt_implicit, args, ,
	   "implicit integration method")
{
#ifdef OCTAVE_DEV
  static 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 ();
  static ColumnVector	dx	= args(1).column_vector_value ();
  static Matrix		AA	= args(2).matrix_value ();
  static 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 // !OCTAVE_DEV
  static 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 ();
  static ColumnVector	dx	= args(1).vector_value ();
  static Matrix		AA	= args(2).matrix_value ();
  static 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 // OCTAVE_DEV
#endif // STANDALONE

  register int row, col;

  for (row = 0; row < Nx; row++)
    {
      if (0 != openx (row))
      if (openx (row) > 0.5)
	{
	  AAx (row) = 0.0;
          dx (row) = 0.0;
#ifndef RISKLEAKS
	  for (col = 0; col < Nx; col++)
	    {
	      AA (row,col) = 0.0;
	      AA (col,row) = 0.0;
	    }
#endif
	}
    }

  x = static_cast<ColumnVector> (xleftdiv (AA, static_cast<Matrix>(AAx + dx * t)));
  
  for (row = 0; row < Nx; row++)
    {
      if (0 != openx (row))
      if (openx (row) > 0.5)
	{
	  x (row) = 0.0;
	}
    }

#ifdef STANDALONE
  return x;