/* Template file for coding d.d.e. models for solv95 - Example is a switched dde
model
*/
#include <math.h>
//#include <windows.h>
#include "solv95.h"
/***************************************************************************/
/* Put global variables here. These should never be written to from */
/* grad() or switchfunctions(), directly or indirectly. */
/***************************************************************************/
/***************************************************************************/
/* Problem specific routines */
/***************************************************************************/
void initcons(no_vars,no_cons) int *no_cons,*no_vars;
/* this routine specifies the number of constants and number of variables in
the model.
*/
{ *no_cons=6;
*no_vars=2;
}
void switchfunctions(sw,s,c,t)
double *sw,*s,*c,t;
/* This routine sets the values of the switch functions. When the switch
functions pass through zero from positive to negative the state variables
may be reset in function map(). The switch functions should pass smoothly
through 0 and should never have both value and first derivative zero. The
same switch must not pass through zero from positive to negative more than
once in an integration timestep. An example of a switch function is:
sw[0]=sin(pi*t/30.0)
which passes through zero every 60 time units. Switches may include state
variables provided the above conditions are met. Note that to get 'Solver'
style switches define twice as many switches and let e.g. sw[1]=-sw[0] */
{ static double pi;
static int first=1;
if (first) {first=0;pi = 2*asin(1.0);}
sw[0]=sin(2*pi*t/c[1]); // add resource
sw[1]=sin(2*pi*(t-c[4])/c[1]); // stop integrator as discontinuity in lagged variable occurs
}
void map(s,c,t,swno)
double *s,*c,t;int swno;
/* This routine is called whenever one of the switch functions passes through
zero. 'swno' is the number of the switch function. The state variables
can be changed discontinuously within this routine. eg:
if (swno==1)
{ s[0]=coeff[1]*(s[0]);}
time and the coefficients should not be changed.
*/
{ if (swno==0)
{ s[0]+=c[2]; // adding resources
} // otherwise nothing required
}
void grad(g,s,c,t,pastvalue)
double *g,*s,*c,t;
double (*pastvalue)();
/* This routine must provide the gradients g for the state variables s.
So ds[i]/dt=g[i]=fi(s,c,t) where c is the coefficient vector. lagged
variables may be accessed here using pastvalue(i,x,j) which returns the
ith (starting at zero) lagged variable at time x, using lag pointer k
(lag pointers are used by pastvalue to store the history buffer location
corresponding to a lag in order to save exectution time. For example if
your code requires lagged varaible 0 at lags T and 2T for each gradient
calculation then it is efficient to obtain these values using:
pastvalue(0,t-T,0) and pastvalue(0,t-2*T,1) rather than
pastvalue(0,t-T,0) and pastvalue(0,t-2*T,0). The latter works, it's just
slower because more time is spent searching for lagged values)
*/
{ g[0]=-c[0]*s[0]*s[1]; // resource
if (t>c[4])
g[1]=c[3]*pastvalue(0,t-c[4],0)*pastvalue(1,t-c[4],0)-c[5]*s[1]; // consumer
else
g[1]=-c[5]*s[1];
}
void storehistory(his,ghis,g,s,c,t)
double *his,*ghis,*g,*s,*c,t;
/* This is the routine in which the values of the history variables at time
t are calculated and put in the array his, along with gradients in ghis,
using state variables s, gradients of s, g, and coefficients c
e.g. if the state variable 2 is history variable 0, you would need the line:
his[0]=s[2];ghis[0]=g[2];
*/
{ his[0]=s[0];ghis[0]=g[0];
his[1]=s[1];ghis[1]=g[1];
}
void statescale(double *scale)
/* In this routine you can set scale factors for error control. For each
state variable the maximum permisable error will be bounded below by the
tolerance multiplied by scale[i]. If you don't supply values then zero will
be used.
Non-zero scale values are useful for variables that start at zero and
leave zero without 3rd order continuity. */
{
}
void initst(s,c,t)
double *s,*c,t;
/* initialise state variables and any global constants here, you can use c */
{ s[0]=0.0;s[1]=1.0;
}
void initout(out) usercontrol *out;
/* This routine is where the output windows are first set up. You have to
fill in the global structure "out", which tells solv95 which
variables to plot and in which windows */
{ int i;
/* how many windows */
out->no_windows=2;
/* how many lines in each window e.g. for two lines in window zero set
out->lines[0]=2 ..... */
out->lines[0]=1;
out->lines[1]=1;
/* which state variables go in which windows and what are they called?
e.g. to put state variable 4 in window 3 as curve 0, and label it `trash':
out->index[4].win=3;out->index[4].cur=0;out->label[4]="trash";
If you don't provide this information for a variable then it isn't plotted
*/
out->label[0]="Resource";
out->index[0].win=0;out->index[0].cur=0;
out->label[1]="Consumer";
out->index[1].win=1;out->index[1].cur=0;
/* labels for windows */
out->wname[0]="Output";
out->wname[1]="More Output";
/* Now set up initial range of y axis variables for each window */
out->range[0].y0=0.0;out->range[0].y1=100.0;
out->range[1].y0=0.0;out->range[1].y1=200.0;
out->xlabel="Time";
/* initialise file output. Set out->fileno to the number of columns to output
(in addition to time). Set out->fout[i] to the number of the state variable
that is to go in column i. Set out->fileno=0 for no output. */
out->fileno=2;
out->fout[0]=0;
out->fout[1]=1;
/* integration details*/
out->t0=0.0; /* default start time */
out->t1=300.0; /* default stop time */
out->dt=1.0; /* initial timestep */
out->tol=0.000005; /* integration tolerance */
out->dout = 1.0; /* approximate average output timestep */
out->hbsize=1000L; /* how many past values to store for each history variable */
out->nhv=2; /* Number of history (lagged) variables */
out->nlag=1; /* Number of lag markers per history variable (set to 1 if unsure)*/
out->nsw=2; /* number of switch varaibles */
/* Initial Values and names for constants to be prompted for at run time */
out->c[0]=0.1; out->cname[0]="alpha";
out->c[1]=10.0; out->cname[1]="T";
out->c[2]=50.0;out->cname[2]="R_a";
out->c[3]=0.05;out->cname[3]="gamma";
out->c[4]=5.0;out->cname[4]="tau";
out->c[5]=0.02;out->cname[5]="delta";
/* Optional information text on parameters - split text between lines with `\'*/
out->cinfo[0]="Consumption rate per consumer per unit resource.";
out->cinfo[1]="Time between resource top up.";
out->cinfo[2]="Amount by which resources are topped up.";
out->cinfo[3]="Consumer per capita fecundity per unit resource";
out->cinfo[4]="Consumer development time.";
out->cinfo[5]="Consumer per capita death rate.";
/* Text and title for initial welcome box (optional) */
out->initialtitle="Switched DDE Model";
out->initialtext=" This model is an example that came with solv95 ";
}
body
{
background-color: #FFF;
font-size: 13px;
font-family: verdana, helvetica, arial, sans-serif;
}
h1, h2, h3, h4, h5, h6
{
color: #000080;
}
typedef struct
{ double y0,y1;} y0y1type;
typedef struct
{ float x,y;} xytype;
typedef struct
{ int win,cur;} wctype;
typedef struct
{ double x0,x1;} x0x1type;
typedef struct
{
int no_windows;
int *lines; /* variables per window */
int nhv;
int nlag;
int nsw;
int fileno;
int *fout;
long hbsize;
wctype *index; /* window and curve for state variable */
char **label;
char **cname;
char **wname;
char *initialtext;
char *initialtitle;
char **cinfo;
char *xlabel;/* labels for state variables and constants*/
y0y1type *range; /* y-range for [window] */
double t0;
double t1;
double dt;
double dout;
double *c;
double tol;
} usercontrol;
void initcons(int * no_vars,int *no_cons);
void initout(usercontrol *out);
void initialstate(double *s, double *c, double t);
Using DDESolve
The main window
Run the programme. The main window will appear to allow you to edit the model parameters:
Double-click on a value to edit: you may not edit the parameter names or information strings.
Click on the run button to run. Any graphs requested by the model code will now be displayed.
The Plot Window Inspector
This can be generated by Graph->Inspector menu item.
Plotting State Variables
The Graph->Plot menu item will bring up the plot panel.
Select whichever state variables you wish to plot by selecting the appropriate rows in the table. Use the pull-down to select the abscissa: t may be time or any of the state variables.
The About Box
This will display the initialtitle and initialtext strings defined by the model. It will also tell you the process ID of the DDESolve application.
Back to index
Bud1 o l v e . a D D E S o l v e . a p pfwi0blob icnv D D E S o l v e . a p pfwswlong D D E S o l v e . a p picgoblob @ @ @ @ E DSDB ` @ @ @ `T T 8__PAGEZERO __TEXT __text __TEXT __symbol_stub __TEXT __picsymbol_stub__TEXT $__symbol_stub1 __TEXT __cstring __TEXT ׀ 't ǀ __literal4 __TEXT < __literal8 __TEXT 0 0 __DATA __data __DATA __nl_symbol_ptr __DATA $ __la_symbol_ptr __DATA h # __dyld __DATA ` ` __cfstring __DATA h h __bss __DATA X __common __DATA __OBJC __cat_cls_meth __OBJC __cat_inst_meth __OBJC __string_object __OBJC __cstring_object__OBJC __message_refs __OBJC __sel_fixup __OBJC __cls_refs __OBJC __class __OBJC __meta_class __OBJC __cls_meth __OBJC __inst_meth __OBJC
__protocol __OBJC # __category __OBJC # __class_vars __OBJC # __instance_vars __OBJC # __module_info __OBJC ) __symbols __OBJC * 8__LINKEDIT 0 P @T /usr/lib/dyld T DvS /System/Library/Frameworks/Cocoa.framework/Versions/A/Cocoa < Dl /usr/local/lib/libplplotd.9.dylib < Dlh /usr/local/lib/libcsirocsa.0.dylib T C, /Library/Frameworks/AquaTerm.framework/Versions/A/AquaTerm 4 DvS /usr/lib/libobjc.A.dylib 8 DvS /usr/X11R6/lib/libX11.6.dylib 4 DvS X /usr/lib/libSystem.B.dylib 8 Dַ /usr/local/lib/libgcc_s.1.dylib 7 ( P 5 5 T B 6 = 4 R 3 B ( |: x8!T! 48 !z 8 ;c W{:|H |<@ != |#x|yx|+x !b <@ $ / }"Kx@ 8= 9)׀H <|i|lxN!H X/ /A 8B H 9" }"Kx | t/ @<@ fx" H 8 / @<@ ;F B ԀB / A |I|LxN!<@ B ܀B / A |I|LxN!H <@ B / A |I|LxN!<` 8 <8cׄH I <| | xN!<` 8 @8c״H )<` 8 88cH <` 8 H8cH <` 8 L8cH <` 8 P8cH A @/ A, 8/ A H/ A L/ A P/ A |I|Lx? ; ;
%@