--- a/main/odepkg/src/odepkg_mexsolver_dopri5.c	Wed Jan 24 20:50:31 2007 +0000
+++ b/main/odepkg/src/odepkg_mexsolver_dopri5.c	Wed Jan 24 20:51:36 2007 +0000
@@ -36,14 +36,18 @@
 #include <mex.h>       /* Needed for all mex definitions */
 #include <string.h>    /* Needed for the memcpy function etc.*/
 #include "odepkgmex.h" /* Needed for the mex extensions */
-#include "odepkgext.h" /* Needed for the odepkg extensions */
+#include "odepkgext.h" /* Needed for the odepkg interfaces */
 #include "dopri5.h"    /* Needed for the solver function */
 
-/* These are the prototypes from this file */
+/* These prototype definitions from functions that are found in this file */
 void fodefun (unsigned n, double x, double *y, double *f);
 void fsolout (long nr, double xold, double x, double* y, unsigned n, int* irtrn);
 void mexFunction (int nlhs, mxArray* plhs[], int nrhs, const mxArray* prhs[]);
 
+/**
+ * TODO
+ * Return value: On success the constant %true, otherwise %false.
+ **/
 void fodefun (unsigned n, double x, double *y, double *f) {
   int  vcnt = 0;
   int  vnum = 0;
@@ -85,26 +89,68 @@
 
 void fsolout (long nr, double xold, double x, double* y, unsigned n, int* irtrn) {
   bool     vsuc = false;
+  unsigned int vref = 0;
+  unsigned int vcnt = 0;
+  unsigned int vnum = 0;
+  double   vdbt = 0.0;
+
+  double  *vdob = NULL;
+  double  *vdbl = NULL;
   mxArray *vtmp = NULL;
   mxArray *vtem = NULL;
   mxArray *vtim = NULL;
+  mxArray *vrtm = NULL;
+  mxArray *vrvl = NULL;
 
   /* Convert the double *y into a selection of variables */
   fy2mxArray (n, y, &vtem);
   vtim = mxCreateDoubleScalar (x);
+
+  /* Save the solution if this is not the initial first call */
   if (nr > 1) fsolstore (1, &vtim, &vtem);
 
-  /* Call plotting function if this is not the initial first call */
+  /* Call output function if this is not the initial first call */
   fodepkgvar (2, "plotfun", &vtmp);
   if (!mxIsEmpty (vtmp) && nr > 1) {
+    fodepkgvar (2, "refine", &vtmp); /* Check for "Refine" option */
+    vdbl = mxGetPr (vtmp);
+    vref = (int)vdbl[0];
+
+    if (vref > 0) { /* If the Refine option is > 0 */
+      vdob = (double *) mxMalloc (n * sizeof (double));
+
+      for (vcnt = 0; vcnt < vref; vcnt++) {
+        vdbt = xold + (double)(vcnt + 1) * ((x - xold) / (vdbl[0] + 2));
+        vrtm = mxCreateDoubleScalar (vdbt);
+        /* Time stamps for approximation: mexPrintf ("%f %f %f\n", xold, vdbt, x); */
+
+        for (vnum = 0; vnum < n; vnum++) vdob[vnum] = contd5(vnum, vdbt);
+        fy2mxArray (n, vdob, &vrvl);
+        /* Time stamps and approx. values: mexPrintf ("%f %f %f\n", vdbt, vdob[0], vdob[1]); */
+        fodepkgplot (vrtm, vrvl, mxCreateString (""));
+      }
+      mxFree (vdob);
+    }
     vsuc = fodepkgplot (vtim, vtem, mxCreateString (""));
-    irtrn[0] = ((int) vsuc) - 1;
+    if (vsuc == false) irtrn[0] = - 1; /* Stop integration ? */
   }
 
-  /* #ifdef __ODEPKGDEBUG__ */
-  /*   mexPrintf ("%ld  %f  %f  %f  %f  %d  %d\n",  */
-  /*     nr, xold, x, y[0], y[1], n, irtrn[0]); */
-  /* #endif */
+  fodepkgvar (2, "eventfun", &vtmp);
+  if (!mxIsEmpty (vtmp)) {
+    fodepkgevent (vtim, vtem, mxCreateString (""), &vtmp);
+    fodepkgvar (3, "EventSolution", NULL);  /* Remove the last events results */
+    fodepkgvar (1, "EventSolution", &vtmp); /* Set the new events results */
+
+    vtem = mxGetCell (vtmp, 0);    /* Check if we have to stop solving */
+    if (mxIsLogicalScalarTrue (vtem)) {
+      fsolstore (3, NULL, NULL);   /* Remove the last solution entry */
+      vtem = mxGetCell (vtmp, 3);  /* Get last row of events solution */
+      vtim = mxCreateDoubleScalar (x);
+      vtmp = mxGetMatrixRow (vtem, mxGetM (vtem) - 1);
+      fsolstore (1, &vtim, &vtmp); /* Set last row of events solution */
+      irtrn[0] = -1;               /* Tell the solver to stop solving */
+    }
+  }
 }
 
 void mexFunction (int nlhs, mxArray* plhs[], int nrhs, const mxArray* prhs[]) { 
@@ -308,12 +354,10 @@
 
   /* Handle the odeoptions structure field: REFINE */
   fodepkgvar (2, "odeoptions", &vtmp);
-  fodepkgvar (2, "defoptions", &vtem);
-  if (!mxIsEqual (mxGetField (vtmp, 0, "Refine"), mxGetField (vtem, 0, "Refine")))
-    mexWarnMsgTxt ("Not yet implemented option \"Refine\" will be ignored");
+  vtem = mxGetField (vtmp, 0, "Refine");
+  fodepkgvar (1, "refine", &vtem);
 #ifdef __ODEPKGDEBUG__
-  else
-    mexPrintf ("ODEPKGDEBUG: The Refine option has not been set\n");
+  mexPrintf ("ODEPKGDEBUG: The Refine option has been set\n");
 #endif
 
   /* Handle the odeoptions structure field: STATS */
@@ -352,12 +396,10 @@
 
   /* Handle the odeoptions structure field: EVENTS */
   fodepkgvar (2, "odeoptions", &vtmp);
-  fodepkgvar (2, "defoptions", &vtem);
-  if (!mxIsEqual (mxGetField (vtmp, 0, "Events"), mxGetField (vtem, 0, "Events")))
-    mexWarnMsgTxt ("Not yet implemented option \"Events\" will be ignored");
+  vtem = mxGetField (vtmp, 0, "Events");
+  fodepkgvar (1, "eventfun", &vtem);
 #ifdef __ODEPKGDEBUG__
-  else
-    mexPrintf ("ODEPKGDEBUG: The Events option has not been set\n");
+  mexPrintf ("ODEPKGDEBUG: The Events option has been set\n");
 #endif
 
   /* Handle the odeoptions structure field: JACOBIAN */
@@ -475,6 +517,20 @@
 /*   fodepkgvar (9, NULL, NULL); */
 /* #endif */
 
+  /* Initialize the function: FODEPKGEVENT */
+/*   vtem = mxCreateCellArray (1, &vone); */
+/*   fodepkgvar (1, "EventSolution", &vtem); */
+  fodepkgvar (2, "eventfun", &vtmp);
+  if (!mxIsEmpty (vtmp)) {
+    vtmp = mxCreateDoubleScalar (vtimeslot[0]);
+    fy2mxArray (vodedimen, vinitvals, &vtem);
+    if (!fodepkgevent (vtmp, vtem, mxCreateString ("init"), NULL)) {
+      fodepkgvar (2, "eventfun", &vtmp);
+      sprintf (vmsg, "Error at initialisation of event function \"%s\"", mxArrayToString (vtmp));
+      mexErrMsgTxt (vmsg);
+    }
+  }
+
   /* Initialize the function: FSOLSTORE */
   vtmp = mxCreateDoubleScalar (vtimeslot[0]);
   fy2mxArray (vodedimen, vinitvals, &vtem);
@@ -569,8 +625,8 @@
       mxAddField (plhs[0], "stats");
       vnum = mxGetFieldNumber (plhs[0], "stats");
       mxSetFieldByNumber (plhs[0], 0, vnum, vtem);
+    }
 
-    }
   }
 
   else if (nlhs == 2) {
@@ -580,7 +636,14 @@
   }
 
   else if (nlhs == 5) {
+    fsolstore (2, &vtmp, &vtem);
+    plhs[0] = vtmp;
+    plhs[1] = mxTransposeMatrix (vtem);
 
+    fodepkgvar (2, "EventSolution", &vtem);
+    plhs[2] = mxDuplicateArray (mxGetCell (vtem, 1));
+    plhs[3] = mxDuplicateArray (mxGetCell (vtem, 2));
+    plhs[4] = mxDuplicateArray (mxGetCell (vtem, 3));
   }
 
 #ifdef __ODEPKGDEBUG__
@@ -596,7 +659,7 @@
   fodepkgvar (2, "plotfun", &vtmp);
   if (!mxIsEmpty (vtmp)) {
     vtmp = mxDuplicateArray (prhs[1]);
-    vtem = mxDuplicateArray (prhs[2]);
+    fy2mxArray (vodedimen, vinitvals, &vtem); /* vtem = mxDuplicateArray (prhs[2]); */
     if (!fodepkgplot (vtmp, vtem, mxCreateString ("done"))) {
       fodepkgvar (2, "plotfun", &vtmp);
       sprintf (vmsg, "Error at finalisation of output function \"%s\"", mxArrayToString (vtmp));
@@ -604,6 +667,18 @@
     }
   }
 
+  /* Initialize the function: FODEPKGEVENT */
+  fodepkgvar (2, "eventfun", &vtmp);
+  if (!mxIsEmpty (vtmp)) {
+    vtmp = mxCreateDoubleScalar (vtimeslot[0]);
+    fy2mxArray (vodedimen, vinitvals, &vtem);
+    if (!fodepkgevent (vtmp, vtem, mxCreateString ("done"), NULL)) {
+      fodepkgvar (2, "eventfun", &vtmp);
+      sprintf (vmsg, "Error at initialisation of event function \"%s\"", mxArrayToString (vtmp));
+      mexErrMsgTxt (vmsg);
+    }
+  }
+
   /* Cleanup all internals of: FSOLSTORE */
   fsolstore (4, NULL, NULL);
   /* Cleanup all internals of: FODEPKGVAR */

--- a/main/odepkg/src/odepkgext.c	Wed Jan 24 20:50:31 2007 +0000
+++ b/main/odepkg/src/odepkgext.c	Wed Jan 24 20:51:36 2007 +0000
@@ -18,23 +18,46 @@
 */
 
 #ifdef HAVE_CONFIG_H
-  #include "config.h"  /* Needed for GCC_ATTR_UNUSED if compiled with mkoctfile */
+  #include "config.h"  /* Needed for #defines if compiled with mkoctfile */
 #endif
 
 #include <mex.h>       /* Needed for all mex definitions */
-#include <string.h>    /* Needed for the strcmp function */
+#include <string.h>    /* Needed for the memcpy function */
 #include <stdio.h>     /* Needed for the sprintf function */
-#include "odepkgext.h"
-#include "odepkgmex.h"
+#include "odepkgext.h" /* Needed for the odepkg interfaces */
+#include "odepkgmex.h" /* Needed for the mex extensions */
 
-bool fodepkgvar (const unsigned int vtodo, const char *vname, mxArray **vvalue) {
+/**
+ * fodepkgvar - Stores mxArrays in an internal variable
+ * @vdeci: The decision for what has to be done
+ * @vname: The name of the variable that is stored
+ * @vvalue: The value of the stored variable
+ *
+ * If @vdeci is 0 (@vname is %NULL and @vvalue is %NULL) then this
+ * function is initialized. Otherwise if @vdeci is 1 (@vname must be
+ * the name of the variable and @vvalue must be a pointer to its
+ * value) the variable is set and will be stored, if @vdeci is 2 the
+ * value of the stored variable will be returned (@vname is the name
+ * of the variable and @vvalue is a pointer to the place to where its
+ * value should be returned). If @vdeci is 3 then the variable is
+ * removed from the internal storage system (@vname must be the name
+ * of the variable and @vvalue is a %NULL pointer). @vdeci should be
+ * set 4 if cleaning up has to be performed normally at the end of a
+ * function call (@vname is %NULL and @vvalue is %NULL).
+ *
+ * This function will be mainly used to suppress the use of global
+ * variables in C-function with mex-support.
+ *
+ * Return value: On success the constant %true, otherwise %false.
+ **/
+bool fodepkgvar (const unsigned int vdeci, const char *vname, mxArray **vvalue) {
   int  cone = 1;
   bool vret = false;
   int  vnmb = 0;
 
   static mxArray *vodesetvar = NULL;
 
-  switch (vtodo) {
+  switch (vdeci) {
     case 0:  /* Initialisation has to be performed */
       vodesetvar = mxCreateStructArray (1, &cone, 0, NULL);
       vret = true;
@@ -64,10 +87,12 @@
       vret = true;
       break;
 
+#ifdef __ODEPKGDEBUG__
     case 9:  /* Disp the whole structure (for debugging purpose) */
       mexCallMATLAB (0, NULL, 1, &vodesetvar, "disp");
       vret = true;
       break;
+#endif
 
     default: /* Nothing has to done, get out of here */
       break;
@@ -75,14 +100,36 @@
   return (vret);
 }
 
-bool fsolstore (unsigned int vz, mxArray **vt, mxArray **vy) {
+/**
+ * fsolstore - Stores a mxArray time- and solution vector
+ * @vdeci: The decision for what has to be done
+ * @vt: A pointer to the time value mxArray
+ * @vy: A pointer to the mxArray solution vector
+ *
+ * If @vdeci is 0 (@vt is a pointer to the initial time step and @vy
+ * is a pointer to the initial values vector) then this function is
+ * initialized. Otherwise if @vdeci is 1 (@vt is a pointer to another
+ * time step and @vy is a pointer to the solution vector) the values
+ * of @vt and @vy are added to the internal variable, if @vdeci is 2
+ * then the newly allocated internal vectors are returned as pointer
+ * @vt and @vy. If @vdeci is 4 (@vdeci being 3 is unused for
+ * compatibility with function fodepkgvar) then the internal variable
+ * is cleaned up and the memory is free.
+ *
+ * This function will be mainly used to suppress the use of global
+ * variables for storing a time and solution vector in C-functions
+ * with mex-support.
+ *
+ * Return value: On success the constant %true, otherwise %false.
+ **/
+bool fsolstore (unsigned int vdeci, mxArray **vt, mxArray **vy) {
   bool vret = false;
   static int vdim = 0;
   static int vcnt = 0;
   static double *vtstore = NULL;
   static double *vystore = NULL;
 
-  switch (vz) {
+  switch (vdeci) {
     case 0:
       vtstore = (double *) mxMalloc (1 * sizeof (double));
       memcpy ((void *) (vtstore + 0), (void *) mxGetPr (vt[0]), 1 * sizeof (double));
@@ -106,7 +153,9 @@
       vret = true; break;
 
     case 3:
-      break;
+      vtstore = (double *) mxRealloc (vtstore, (vcnt-1) * sizeof (double));
+      vystore = (double *) mxRealloc (vystore, (vcnt-1) * vdim * sizeof (double));
+      vcnt--; vret = true; break;
 
     case 4:
       mxFree (vtstore);
@@ -120,6 +169,21 @@
   return (vret);
 }
 
+/**
+ * fy2mxArray - Takes elements of a vector and returns a mxArray
+ * @n: The number of elements in the vector @y
+ * @y: A pointer to the first element of a vector
+ * @vval: A pointer to the mxArray solution vector
+ *
+ * Takes elements out of a vector @y that has the size @n. The
+ * elements that have to be taken are stored via the fodepkgvar()
+ * function in the variable @outputsel (make sure calling
+ * fodepkgvar(1, "outputsel", @value) before calling this
+ * fy2mxArray()).  The solution is stored in a newly allocated mxArray
+ * and is returned as a pointer @val;
+ *
+ * Return value: On success the constant %true, otherwise %false.
+ **/
 bool fy2mxArray (unsigned int n, double *y, mxArray **vval) {
   int      vnum = 0;
   int      vcnt = 0;
@@ -146,7 +210,25 @@
   return (true);
 }
 
-bool fodepkgplot (mxArray *vtime, mxArray *vvalues, mxArray *vflag) {
+/**
+ * fodepkgplot - Interface to the output function
+ * @vtime: The time value mxArray
+ * @vvalues: The solution vector mxArray
+ * @vdeci: The decision for what has to be done
+ *
+ * Prepares the right-hand side variable before calling the odepkg
+ * output function. @vdeci can either be "init", "" or "done
+ * (cf. interpreter help text eg. odeprint(), odeplot), @vtime is the
+ * actual timevalue and vvalues the solution vector. The output
+ * function that has to be called from fodepkgplot() is set via the
+ * command odeset() at interpreter side. Before calling this function
+ * fodepkgvar(1, "varargin", value) and fodepkgvar(1, "plotfun",
+ * value) have to be called to set the necessary name of the output
+ * function and further arguments that have to be passed.
+ *
+ * Return value: %false if solving has to be stopped, otherwise %true.
+ **/
+bool fodepkgplot (mxArray *vtime, mxArray *vvalues, mxArray *vdeci) {
   bool vret = false;
   int  vcnt = 0;
   int  velm = 0;
@@ -169,7 +251,7 @@
   vrhs[1] = vvalues;
 
   /* Set the vrhs[2] element for the plotting function */
-  if (!mxIsEmpty (vflag)) vrhs[2] = vflag;
+  if (!mxIsEmpty (vdeci)) vrhs[2] = vdeci;
   else vrhs[2] = mxCreateString ("");
 
   /* Set the vrhs[3..] element for the plotting function */
@@ -180,8 +262,9 @@
 
   /* Call the plotting function */
   fodepkgvar (2, "plotfun", &vtmp);
-  if (strcmp (mxArrayToString (vflag), "init") ||
-      strcmp (mxArrayToString (vflag), "done")) {
+  if ((strcmp (mxArrayToString (vdeci), "init") == 0) ||
+      (strcmp (mxArrayToString (vdeci), "done") == 0)) {
+    /* mexPrintf ("-----> INIT or DONE\n"); */
     if (mexCallMATLAB (0, NULL, 3+velm, vrhs, mxArrayToString (vtmp))) {
       sprintf (vmsg, "Calling \"%s\" has failed", mxArrayToString (vtmp));
       mexErrMsgTxt (vmsg);
@@ -189,11 +272,12 @@
     vret = true;
   }
   else {
+    /* mexPrintf ("-----> CALC\n"); */
     if (mexCallMATLAB (1, vlhs, 3+velm, vrhs, mxArrayToString (vtmp))) {
       sprintf (vmsg, "Calling \"%s\" has failed", mxArrayToString (vtmp));
       mexErrMsgTxt (vmsg);
     }
-    if (!mxIsLogicalScalarTrue (vlhs[0])) vret = true;
+    if (mxIsLogicalScalarTrue (vlhs[0])) vret = true;
     else vret = false;
   }
 
@@ -202,6 +286,88 @@
   return (vret);
 }
 
+/**
+ * fodepkgevent - Interface to the event function
+ * @vtime: The time value mxArray
+ * @vvalues: The solution vector mxArray
+ * @vdeci: The decision for what has to be done
+ * @vval: The event output values
+ *
+ * Prepares the right-hand side variable before calling the odepkg
+ * event function. @vdeci can either be "init", "" or "done
+ * (cf. interpreter help text eg. odeprint(), odeplot), @vtime is the
+ * actual timevalue and vvalues the solution vector. The event
+ * function that has to be called from fodepkgplot() is set via the
+ * command odeset() at interpreter side. Before calling this function
+ * fodepkgvar(1, "varargin", value) and fodepkgvar(1, "eventfun",
+ * value) have to be called to set the necessary name of the output
+ * function and further arguments that have to be passed.
+ *
+ * Return value: %false if solving has to be stopped, otherwise %true.
+ **/
+bool fodepkgevent (mxArray *vtime, mxArray *vvalues, mxArray *vdeci, mxArray **vval) {
+  const char *odepkg_event_handle = "odepkg_event_handle";
+
+  bool vret = false;
+  int  vcnt = 0;
+  int  velm = 0;
+  char vmsg[64] = "";
+
+  mxArray  *vtmp = NULL;
+  mxArray **vlhs = NULL;
+  mxArray **vrhs = NULL;
+
+  /* Get number of varargin elements for allocating memory */
+  fodepkgvar (2, "varargin", &vtmp);
+  velm = mxGetNumberOfElements (vtmp);
+  vlhs = (mxArray **) mxMalloc (sizeof (mxArray *));
+  vrhs = (mxArray **) mxMalloc ((4 + velm) * sizeof (mxArray *));
+
+  /* Set the vrhs[0] element for the event function */
+  fodepkgvar (2, "eventfun", &vtmp);
+  vrhs[0] = vtmp;
+
+  /* Set the vrhs[1] element for the event function */
+  vrhs[1] = vtime;
+
+  /* Set the vrhs[2] element for the event function */
+  vrhs[2] = vvalues;
+
+  /* Set the vrhs[3] element for the event function */
+  if (!mxIsEmpty (vdeci)) vrhs[3] = vdeci;
+  else vrhs[3] = mxCreateString ("");
+
+  /* Set the vrhs[4..] element for the event function */
+  fodepkgvar (2, "varargin", &vtmp);
+  if (velm > 0) /* Fill up vrhs[3..] */
+    for (vcnt = 0; vcnt < velm; vcnt++)
+      vrhs[4+vcnt] = mxGetCell (vtmp, vcnt);
+
+  /* Call the event function */
+  if ((strcmp (mxArrayToString (vdeci), "init") == 0) ||
+      (strcmp (mxArrayToString (vdeci), "done") == 0)) {
+    //mexPrintf ("---> BIN ICH DRINNN? %s %d %d\n", mxArrayToString (vdeci), !(-1), !(-2));
+    if (mexCallMATLAB (0, NULL, 4+velm, vrhs, odepkg_event_handle)) {
+      sprintf (vmsg, "Calling \"%s\" has failed", odepkg_event_handle);
+      mexErrMsgTxt (vmsg);
+    }
+    vret = true;
+  }
+  else {
+    //mexPrintf ("---> BIN WOANDERS? %s %d %d\n", mxArrayToString (vdeci), !(-1), !(-2));
+    if (mexCallMATLAB (1, vlhs, 4+velm, vrhs, odepkg_event_handle)) {
+      sprintf (vmsg, "Calling \"%s\" has failed", odepkg_event_handle);
+      mexErrMsgTxt (vmsg);
+    }
+    for (vcnt = 0; vcnt < 3; vcnt++) vval[vcnt] = vlhs[vcnt];
+    vret = true;
+  }
+
+  mxFree (vlhs); /* Free the vlhs mxArray */
+  mxFree (vrhs); /* Free the vrhs mxArray */
+  return (vret);
+}
+
 /*
 Local Variables: ***
 mode: C ***

--- a/main/odepkg/src/odepkgext.h	Wed Jan 24 20:50:31 2007 +0000
+++ b/main/odepkg/src/odepkgext.h	Wed Jan 24 20:51:36 2007 +0000
@@ -20,11 +20,11 @@
 #ifndef __ODEPKGEXT__
 #define __ODEPKGEXT__ 1
 
-extern bool fodepkgvar (const unsigned int vtodo, const char *vname, mxArray **vvalue);
-extern bool fsolstore (unsigned int vz, mxArray **vt, mxArray **vy);
-extern bool fy2mxArray (unsigned int n, double *y, mxArray **vval);
-extern bool fodepkgplot (mxArray *vtime, mxArray *vvalues, mxArray *vflag);
-
+extern bool fodepkgvar   (const unsigned int vdeci, const char *vname, mxArray **vvalue);
+extern bool fsolstore    (unsigned int vdeci, mxArray **vt, mxArray **vy);
+extern bool fy2mxArray   (unsigned int n, double *y, mxArray **vval);
+extern bool fodepkgplot  (mxArray *vtime, mxArray *vvalues, mxArray *vdeci);
+extern bool fodepkgevent (mxArray *vtime, mxArray *vvalues, mxArray *vdeci, mxArray **vval);
 #endif /* __ODEPKGEXT__ */
 
 /*