How to update an old mex code to work again with ode45?
1 view (last 30 days)
Show older comments
Hello
Back in 94 the following mex code worked with ode45 (Matlab 4.2 on solaris 2.6). Although I have the feeling that this type of code has been completely deprecated several years ago, I wonder whether it is possible to translate it to a new version (matlab 2017a). (I know that a m-file would do the job but I have lots of similar codes ..).
/* Lorenz Equation (84) - Simulation */
#include "mex.h"
/* #include "mymex.h" */
/* Define the system sizes here: */
#define NSTATES 12
#define NOUTPUTS 12
#define NCOEFFS 16
static Matrix *Coeffs[NCOEFFS];
#define a mxGetPr(Coeffs[0])
#define b mxGetPr(Coeffs[1])
#define F mxGetPr(Coeffs[2])
#define G mxGetPr(Coeffs[3])
#define X0 mxGetPr(Coeffs[4])
#define X1 mxGetPr(Coeffs[5])
#define X2 mxGetPr(Coeffs[6])
#define X01 mxGetPr(Coeffs[7])
#define X11 mxGetPr(Coeffs[8])
#define X21 mxGetPr(Coeffs[9])
#define X02 mxGetPr(Coeffs[10])
#define X12 mxGetPr(Coeffs[11])
#define X22 mxGetPr(Coeffs[12])
#define X03 mxGetPr(Coeffs[13])
#define X13 mxGetPr(Coeffs[14])
#define X23 mxGetPr(Coeffs[15])
static
#ifdef __STDC__
void init_conditions(double *x0)
#else
void init_conditions(x0)
double *x0;
#endif /* __STDC__ */
{
x0[0] = X0[0];
x0[1] = X1[0];
x0[2] = X2[0];
x0[3] = X01[0];
x0[4] = X11[0];
x0[5] = X21[0];
x0[6] = X02[0];
x0[7] = X12[0];
x0[8] = X22[0];
x0[9] = X03[0];
x0[10] = X13[0];
x0[11] = X23[0];
}
static
#ifdef __STDC__
void derivatives(double t, double *x, double *u, double *dx)
#else
void derivatives(t,x,u,dx)
double t, *x, *u; /* Input variables */
double *dx; /* Output variable*/
#endif /* __STDC__ */
{
int i;
dx[0] = -x[1] * x[1] - x[2] * x[2] - a[0] * (x[0] - F[0]);
dx[1] = x[0] * x[1] - b[0] * x[0] * x[2] - x[1] - G[0];
dx[2] = b[0] * x[0] * x[1] + x[0] * x[2] - x[2];
for (i=0;i<=2;i++)
{
dx[3+i] = (-a[0]) * x[3+i] + (-2.0 * x[1]) * x[6+i] + (-2.0 * x[2]) * x[9+i];
dx[6+i] = (x[1] - b[0] * x[2]) * x[3+i] + (x[0] - 1.0) * x[6+i] + (-b[0] * x[0]) * x[9+i];
dx[9+i] = (b[0] * x[1] + x[2]) * x[3+i] + (b[0] * x[0]) * x[6+i] + (x[0] - 1.0) * x[9+i];
}
}
static
#ifdef __STDC__
void outputs(double t, double *x, double *u, double *y)
#else
void outputs(t,x,u,y)
double t, *x, *u; /* Input variables */
double *y; /* Output variable*/
#endif /* __STDC__ */
{
int i;
for (i=0;i<NOUTPUTS;y[i]=x[i],i++);
}
/* #include "simulink.h" */
The function was called by issuing
[t,x,y]=ode45('lor84',tf,[],options,[],a,b,F,G,X0,X1,X2,X3,X4,X5,X6,X7,X8,X9,X10,X11);
Many thanks.
Ed
4 Comments
Answers (0)
See Also
Categories
Find more on Environment and Settings in Help Center and File Exchange
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!