I'm trying to recreate computation of a SIR model as described here, with extra midpoint calculations. But for some reason no values actually change during the Euler calculations.
#include <iostream>
using namespace std;
int main(){
double s[100];
double i[100];
double r[100];
double bb = 1/2;
double kk = 1/3;
s[0] = 1;
i[0] = 1.27e-6;
r[0] = 0;
int h = 1;
int max = 20;
for (int j = 0; j < max; j += h){
s[j + 1] = s[j] - h*(bb * s[j] * i[j]);
i[j + 1] = i[j] + h*(bb * s[j] * i[j] - kk * i[j]);
r[j + 1] = r[j] + h*(kk * i[j]);
cout << j << "\t" << s[j] << "\t" << i[j] << "\t" << r[j] << "\n";
/*
s[j + 1] = s[j] - 0.5*h*(bb * s[j] * i[j] + bb * s[j+1] * i[j+1]);
i[j + 1] = i[j] + 0.5*h*(bb * s[j] * i[j] - kk * i[j] + bb * s[j+1] * i[j+1] - kk * i[j+1]);
r[j + 1] = r[j] + 0.5*h*(kk * i[j] + kk * i[j+1]);
*/
}
}
Your variables bb
and kk
are both zero due to integer division. Always use double literals:
double bb = 1.0/2.0;
double kk = 1.0/3.0;