/* Copyright (C) 1991-2000 Simon N. Wood  snw@st-and.ac.uk

This program 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
of the License, or (at your option) any later version.

This program 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.
(www.gnu.org/copyleft/gpl.html)

You should have received a copy of the GNU General Public License
along with this program; if not, write to the Free Software
Foundation, Inc., 59 Temple Place - Suite 330, Boston, MA  02111-1307,
USA.*/

/**************************************************************************/
/* Routines to allow DDEfit to more efficiently fit discrete time models. */
/* At run time DDEfit establishes (according to user supplied information */
/* from the model template) whether the model is disctrete time or not.   */
/* If it is then dde() is set to ddem().                                  */
/**************************************************************************/



#include <math.h>
#include <stdlib.h>
#include <stdio.h>
#include "spline.h"
#include "ddeq.h"

void ddem(s,c,t0,t1,dt,eps,dout,ns,nsw,nhv,hbsize,nlag,reset,step)
double *s,      /* State variables */
       *c,      /* coefficients */
       t0,t1,   /* start and stop times */
       *dt,     /* pointer to initial timestep (returns final step - which
		   is step that would have been used if t1 not reached!) */
       eps,     /* fractional tolerance for adaptive stepping */
       dout;    /* interval for output via user routine output(). Every
		   time a step passes 1 or more times of the form t0+i*dout
		   output() is called once. Hence output is only roughly
		   regular. dout=0.0 for no output. */
long hbsize;    /* The number of elements to store in the history buffer */
int nsw,        /* number of switches */
    ns,         /* number of state variables */
    nhv,        /* number of lagged variables */
    reset,      /* set to 0 not to reset, to 1 to reset */
    nlag,       /* number of place markers per history variable */
    step;       /* 0 - adaptive stepping; 1 - record adaptive; 2 - playback
		             adaptive; 3 - free record memory and use adaptive;
		             4 - fixed step */

/* this routine steps a difference equation.... hence the only arguments that
   matter are s,c,t0,t1 - the others simply maintain consistency with
   dde solving code. Input state is state at t0. */


{ double t;
  t=(double)round(t0);
  if (reset&&(dout>0.0))
  poutput(s,t,(void *)NULL,0);
  while (t<t1)
  { map(s,c,t,0); /* t -> t+1 */
    t+=1.0;t=(double)round(t);
    if (dout>0.0) poutput(s,t,(void *)NULL,0);
  }
}


int is_model_discrete(int i)
// routine to signal that a discrete time model is being used here
// call with 0 to query, 1 to set for discrete model, 2 to set for continuous model 
{ static j=0;  
  if (!i) return(j);
  else if (i==1) j=1; else if (i==2) j=0;
  return(j);
}