#! /bin/sh

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

###############################################################
## Version control history
###############################################################
## $Id$
## $Log$




## Revision 1.74  2002/05/01 17:30:56  geraint
## Improved pre-processor directives to better accommodate future alternatives (matlab)
## if necessary.
##
## Revision 1.73  2002/05/01 12:24:41  geraint
## Removed unnecessary inclusion of load-save.h.
##

#include "${sys}_def.h"
#include "${sys}_sympar.h"

#include "mtt_${algebraic_solver}.hh"

#include <fstream>
#if (CODEGENTARGET == STANDALONE)
#include <csignal>

extern ColumnVector ${sys}_ae (
	ColumnVector &x,
	ColumnVector &u,
	const double &t,
	ColumnVector &par);

	const ColumnVector &open_switches);
 
EOF
    ;;
esac
cat <<EOF >> $filename

void set_signal_handlers (void);

#endif // (CODEGENTARGET == STANDALONE)

ColumnVector
mtt_ae (ColumnVector &x,
	ColumnVector &u,
	const double &t,
	ColumnVector &par)

  set_global_value ("MTT_data", data);
#endif // (CODEGENTARGET == STANDALONE)
    row = 0;
    fcputime.close();
    }
}

#if (CODEGENTARGET == STANDALONE)
void dump_data (std::ostream &file)
{
  ColumnVector null (0);
  mtt_write (0.0, null, null, 0, true, file);
}

void handle_signal (int signum)
{
  // handle some signals to ensure data is written.
  std::cerr << "# Writing data to MTT.core (signal " << signum << ")" << std::endl;
  std::ofstream corefile ("MTT.core");
  dump_data (corefile);
  switch (signum)
    {
      case SIGFPE:
        // Intel chips do not raise SIGFPE for DIVZERO :-(
//        raise (SIGABRT);
        break;
      case SIGINT:
	break;
      case SIGQUIT:
	signal (SIGQUIT, SIG_DFL);
	raise (SIGQUIT);
	break;
      default:
        std::cerr << "# Warning: make_ode2odes needs updating!" << std::endl;
	signal (signum, SIG_DFL);
	raise (signum);
        break;
    }
  corefile.close ();
  set_signal_handlers ();
}

void set_signal_handlers (void)
{
  signal (SIGFPE,	handle_signal);
  signal (SIGINT,	handle_signal);
  signal (SIGQUIT,	handle_signal);
}

int main (void) {
  set_signal_handlers ();
#elif (CODEGENTARGET == OCTAVEDLD)
DEFUN_DLD (${sys}_ode2odes, args, ,
"Octave ode2odes representation of system with $method integration method\nUsage: mtt_data = ${sys}_ode2odes (x0, par, simpar)\n")
{
  static octave_value_list retval;
#endif // (CODEGENTARGET == STANDALONE)
  static ColumnVector	x0;
  static ColumnVector	par;
  static ColumnVector	simpar;

  static double
    first	= 0.0,
    dt		= 0.0,
    last	= 0.0,
    stepfactor	= 0.0;
#if (CODEGENTARGET == OCTAVEDLD)
  int nargin = args.length ();
  switch (nargin)
    {
    case 3:
      first		= args (2).${vector_value} ()(0);
      dt		= args (2).${vector_value} ()(1);
      last		= args (2).${vector_value} ()(2);
      stepfactor     	= args (2).${vector_value} ()(3);
      par    		= args (1).${vector_value} ();
      x0      		= args (0).${vector_value} ();
      break;
    case 2:
      first		= mtt_simpar ()(0);
      dt		= mtt_simpar ()(1);
      last		= mtt_simpar ()(2);
      stepfactor     	= mtt_simpar ()(3);
      par    		= args (1).${vector_value} ();
      x0      		= args (0).${vector_value} ();
      break;
    case 1:
      first		= mtt_simpar ()(0);
      dt		= mtt_simpar ()(1);
      last		= mtt_simpar ()(2);
      stepfactor     	= mtt_simpar ()(3);
      par    		= mtt_numpar ();
      x0      		= args (0).${vector_value} ();
      break;
    case 0:
#endif // (CODEGENTARGET == OCTAVEDLD)
      first		= mtt_simpar ()(0);
      dt		= mtt_simpar ()(1);
      last		= mtt_simpar ()(2);
      stepfactor     	= mtt_simpar ()(3);
      par    		= mtt_numpar ();
      x0      		= mtt_state (par);
#if (CODEGENTARGET == OCTAVEDLD)
      break;
    default:
      usage("${sys}_ode2odes (x par simpar)", nargin);
      error("aborting.");
    }
#endif // (CODEGENTARGET == OCTAVEDLD)
  static ColumnVector	dx (MTTNX, 0.0);
  static ColumnVector	x (MTTNX, 0.0);
  static ColumnVector	u (MTTNU, 0.0);
  static ColumnVector	y (MTTNY, 0.0);

  static Matrix	AA (MTTNX, MTTNX, 0.0);
  static ColumnVector	AAx (MTTNX, 0.0);

  for (register int i = 0; i < MTTNY; i++)
    {
      y (i) = 0.0;
    }

  for (register int i = 0; i < MTTNX; i++)
    {
      x (i) = x0 (i);
    }

  for (register int j = 0, i = 1; i <= ilast; i++)
    {
      u	= mtt_input (x, y, t, par);
      y	= mtt_output (x, u, t, par);
      if (0 == j)
	{
           mtt_write (t, x, y, nrows);
	}
EOF
case "$method" in
    "rk4")

          t2 = t + ddt;

        ColumnVector
          x1 (x),
          x2 (x),
          x3 (x);

        k1 = ddt * mtt_rate (x , u, t , par); x1 += k1 * 0.5;
        k2 = ddt * mtt_rate (x1, u, t1, par); x2 += k2 * 0.5;
        k3 = ddt * mtt_rate (x2, u, t1, par); x3 += k3;
        k4 = ddt * mtt_rate (x3, u, t2, par);
        dx = (k1 + 2.0 * (k2 + k3) + k4) / (6.0 * ddt);
      }
EOF
    ;;
    "dassl")
    ;;
    "implicit")
	cat << EOF >> $filename
      dx = mtt_rate (x, u, t, par);
      AA = mtt_smxa (x, u, ddt, par);
      AAx = mtt_smxax (x, u, ddt, par);
EOF
    ;;
    "euler" | *)
	cat << EOF >> $filename
      dx = mtt_rate (x, u, t, par);
EOF
    ;;
esac

## Common stuff
cat <<EOF >> $filename
      open_switches = mtt_logic (x, u, t, par);
      x =  $algorithm; 
      t += ddt;
      j++;
      j = (j == static_cast<int> (stepfactor)) ? 0 : j;
    }



#if (CODEGENTARGET == STANDALONE)


























































  return 0;

#elif (CODEGENTARGET == OCTAVEDLD)




/*  



  retval (0) = octave_value (y);




  retval (1) = octave_value (x);






















  retval (2) = octave_value (t);

*/
  retval = octave_value (get_global_value ("MTT_data"));
  return (retval);

#endif // (CODEGENTARGET == STANDALONE)
}
EOF
}

case ${lang} in
    m)
        make_m	
	;;
    cc)
	make_cc
	;;
    *)
	echo Language ${lang} is not supported
esac