Public feed
#726483
#include <stdlib.h>
#include <stdio.h>
#include "math.h"
#include <mex.h>
#define n 20.0f
#define n_ 3.0f
#define mu 0.10f
#define mu_ 14.0f
#define gamma 0.10f
#define sigma 0.01f
#define beta 0.50f
#define beta_ 0.50f
#define STEPS 0.01f
#define TEND 365.0f
#define COUNT (int)(TEND / STEPS)
int i;
double yvec[COUNT][2];
void getResults()
{
yvec[0][0] = 0.0f; yvec[0][1] = 0.0f;
#define deltay (double)(sigma*(mu*n + gamma*yvec[i][0]) + beta*(n-yvec[i][0])*(yvec[i][1]/n_) - (mu + gamma)*yvec[i][0])
#define deltay_ (double)(beta_*yvec[i][0]*(1.0f - (yvec[i][1]/n_)) - mu_*yvec[i][1])
for (i = 1;i < COUNT; i++) // Euler's Method
{
yvec[i+1][0] = yvec[i][0] + STEPS * deltay;
yvec[i+1][1] = yvec[i][1] + STEPS * deltay_;
}
return;
}
double * myvec;
void mexFunction(int nlhs, mxArray *plhs[],
int nrhs, const mxArray *prhs[])
{
getResults();
plhs[0] = mxCreateDoubleMatrix(COUNT, 2, mxREAL);
myvec = mxGetPr(plhs[0]);
myvec = yvec; //this line is wrong, what should it really be?
return;
}