Plot in Improved Euler's Method problem [closed]

Question

It is given the following problem: $$x'=rx-bx^2-cxy-\frac{\beta xy}{\alpha+x}$$ $$y'=-\mu y+\frac{\beta xy}{\alpha+x}$$ With parameters $r=1,b=1,c=0.01, \mu=0.4,\alpha= 0.36$ and initial values $x(0)=0.2, y(0)=0.05$

I have solved this problem with the NDSolve command

r = 1;
b = 1;
c = 0.01;
a = 0.36;
\[Mu] = 0.4;
deq1 = x'[t] == 
  r* x[t] - b*x[t]^2 - c*x[t]*y[t] - 0.75*x[t]*y[t]/(a + x[t]) 
deq2 = y'[t] == -\[Mu]*y[t] + 0.75*x[t]*y[t]/(a + x[t]) 
soln = NDSolve[{deq1, deq2,  x[0] == 0.2, y[0] == 0.05}, {x[t], y[t]}, {t, 0, 200}]
soln[t_] = 
   NDSolveValue[{deq1, deq2, x[0] == 0.2, y[0] == 0.05}, {x[t], 
       y[t]}, {t, 0, 200}];
Plot[soln[t], {t, 0, 200}, PlotRange -> All]

Now, I am trying the same with the Improved Euler's Method

UPDATED

r = 1;
b = 1;
c = 0.01;
a = 0.36;
\[Mu] = 0.4;
f[x_,y_]:=r* x[t] - b*x[t]^2 - c*x[t]*y[t] - 0.75*x[t]*y[t]/(a + x[t])
g[x_,y_]:=y'[t] == -\[Mu]*y[t] + 0.75*x[t]*y[t]/(a + x[t]) 
Q[a_, b_, h_, N_] := (u[0] = a; v[0] = b; 
  Do[{u[n + 1] = 
     u[n] + h*
       F[u[n] + (h/2)*F[u[n], v[n]], 
        v[n] + (h/2)*
          G[u[n], v[
            n]]],                                                     \
          
    v[n + 1] = 
     v[n] + h*
       G[u[n] + (h/2)*F[u[n], v[n]], v[n] + (h/2)*G[u[n], v[n]]]}, {n,
     0, N}])
Q[0.2, 0.05, 1,200]
u[200]
X = Interpolation[Table[{0.05 n, u[n]}, {n, 0, 100}]]
Y = Interpolation[Table[{0.05 n, v[n]}, {n, 0, 100}]]
ParametricPlot[{X[t], Y[t]}, {t, 0, 200}]

I have run this code on Mathematica but it did not give me results for the plot.

Answer 1

There are several small issues in your code. I would list the most important:

• in your g[x_, y_] definition you forgot to delete y'[t], (as you properly did for f)
• in the rhs of your f and g, x and y are the variables -- in actuality they will be elements of a vector, so the dependence on t should be erased: use just x and y not x[t], y[t].
• your Euler formula uses F and G, while the functions you defined are f and g
• there are some inconsistencies with regards to the maximum t for your integration.

The code below, based on yours, attempts (and succeeds) to create the same graphic as the NDSolve version:

r = 1;
b = 1;
c = 0.01;
a = 0.36;
\[Mu] = 0.4;
f[x_,y_]:=r* x - b*x^2 - c*x*y - 0.75*x*y/(a + x);
g[x_,y_]:= -\[Mu]*y + 0.75*x*y/(a + x);
Q[a_, b_, h_, N_] := (
  u[0] = a; v[0] = b; 
    Do[{
       u[n + 1] =  u[n] + h *f[u[n] + h/2*f[u[n], v[n]],v[n] + h/2* g[u[n], v[n]]],                                                       
       v[n + 1] =   v[n] + h  *g[u[n] + h/2*f[u[n], v[n]], v[n] + h/2*g[u[n], v[n]]]
    }, {n,0, N-1}]);
Q[0.2, 0.05, 1,200];
X = Interpolation[Table[{ n, u[n]}, {n, 0, 200}]];
Y = Interpolation[Table[{n, v[n]}, {n, 0, 200}]];
Plot[{X[t], Y[t]}, {t, 0, 200}]