MATLAB implementation-Prey Model

 

Repeat the problem:

An example of the nonlinear differential equation of the Prey-to-Prey Model is the prey-to-Prey Model. Let X (t) and Y (t) represent the numbers of rabbits and foxes at the moment t, respectively. The Prey-by-prey model shows that x (t) and Y (t) meet the following requirements:

A typical computer simulation coefficient of action:

A = 2, B = 0.02, c = 0.0002, D = 0.8

And meet the following requirements:

X (0) = 3000Rabbit,Y (0) = 120Fox.

In the interval[0, 5]Internal:

(1) UseEulerMethod,HeunMethod and classicLevel 4 Runge-KuttaMethodMultiple step sizes.

SetX=X(T)AndY=Y(T)Draw and explain on the same graph.

(2) For the test variableXAndYThe relationship between them.X,YCoordinatesTo explain the results.

 

Principles:

Euler method

For

DX = f (t, x, y) dt and dy = g (t, x, y) dt

The above formula is used to get:

Divide the interval into M subintervals, H = (B-a)/m, and the grid point is

Then we get the Euler method:

Heun Method

On the interval [,], the points of the differential equation:

Use Trapezoidal Rule to get the following information:

But the right side contains,Use Euler method to obtain heun method:

In each step, the Euler method is used as the prediction, and the Trapezoidal Rule is corrected to obtain the final value.Heun MethodThe general steps are as follows:

Runge-Kutta method of order 4

Based on the derivation process of heum method, we obtain

To improve the accuracy, get the 4-price Runge-Kutta method:

Where:

The following is a program: function E = Euler (F, G, A, B, XA, ya, m) % input-F is the function entered as a string 'F' %-G is the function entered as a string 'G' %-A and B are the left and right endpoints %- XA is the initial condition X () %-Ya is the initial condition y () %-m is the number of steps % output-E = [t 'X'] where T is the vector of abscissas % and X and Y are the vectors of ordinatesh = (b-) /m; t = zeros (1, m + 1); X = zeros (1, m + 1); y = zeros (1, m + 1); t =: h: B; X (1) = XA; y (1) = ya; for j = 1: MX (J + 1) = x (j) + H * feval (F, T (J), x (J), y (j); y (J + 1) = y (j) + H * feval (G, T (J), x (J), y (j); Ende = [t 'X'];

 

function H=heun(f,g,a,b,xa,ya,M)%Input    - f is the function entered as a string 'f'%         - g is the function entered as a string 'g'%             - a and b are the left and right endpoints%             - xa is the initial condition x(a)%             - ya is the initial condition y(a)%             - M is the number of steps%Output - H=[T' X' Y'] where T is the vector of abscissas%                and X and Y are the vectors of ordinatesh=(b-a)/M;T=zeros(1,M+1);X=zeros(1,M+1);Y=zeros(1,M+1);T=a:h:b;X(1)=xa;Y(1)=ya;for j=1:M    k1=feval(f,T(j),X(j),Y(j));    k2=feval(g,T(j),X(j),Y(j));    k3=feval(f,T(j+1),X(j)+h*k1,Y(j)+h*k2);    k4=feval(g,T(j+1),X(j)+h*k1,Y(j)+h*k2);    X(j+1)=X(j)+(h/2)*(k1+k3);    Y(j+1)=Y(j)+(h/2)*(k2+k4);endH=[T' X' Y'];

 

function R=rk4(f,g,a,b,xa,ya,M)%Input    - f is the function entered as a string 'f'%         - g is the function entered as a string 'g'%             - a and b are the left and right endpoints%             - xa is the initial condition x(a)%             - ya is the initial condition y(a)%             - M is the number of steps%Output -R=[T' X' Y'] where T is the vector of abscissas%              and X and Y are the vectors of ordinatesh=(b-a)/M;T=zeros(1,M+1);X=zeros(1,M+1);Y=zeros(1,M+1);T=a:h:b;X(1)=xa;Y(1)=ya;for j=1:M    k1=h*feval(f,T(j),X(j),Y(j));    k2=h*feval(g,T(j),X(j),Y(j));    k3=h*feval(f,T(j)+h/2,X(j)+k1/2,Y(j)+k2/2);    k4=h*feval(g,T(j)+h/2,X(j)+k1/2,Y(j)+k2/2);    k5=h*feval(f,T(j)+h/2,X(j)+k3/2,Y(j)+k4/2);    k6=h*feval(g,T(j)+h/2,X(j)+k3/2,Y(j)+k4/2);     k7=h*feval(f,T(j)+h,X(j)+k5,Y(j)+k6);    k8=h*feval(g,T(j)+h,X(j)+k5,Y(j)+k6);    X(j+1)=X(j)+(k1+2*k3+2*k5+k7)/6;    Y(j+1)=Y(j)+(k2+2*k4+2*k6+k8)/6;endR=[T' X' Y'];

The previous principles have been quite clear, so the computing in the program will not be repeated.

 

Below is the M file:

 

Partial explanation:

 

 

In the same graph, the following is the Y-x image corresponding to Euler, heun, and Runge-Kutta method:

 

Since we can see that the two represent a closed curve in the same figure, it illustrates the self-regulation process of the ecosystem and reveals the dynamic relationship between the number of predators and the number of prey, the close of the image also demonstrates the correctness of the model: the number of images in the cycle is cyclical, and the number of images in the cycle shows an oscillating relationship. In a period of time, this shift is long and increases together for a period of time, however, the growth rate decreases as the number of Parties increases. Specifically, the critical point in both the prey and the prey model is (D/C, A/B) = (4000,100), and it is known that this critical point is the center. For more information, see. The closer the Initial Value Condition x0 is (4000,100), the more the periodic solution is like the elliptic solution of the corresponding linear equations. The feature value of G' (4000,100) is λ = ± √ (AD) * I = ± 2/5 * √ 10i, so the period is (√ 10) π/2, or approximately 4.97

 

Summary:

We can see that in this model, the heun method and Runge-Kutta method have a higher degree of fitting and accuracy than the Euler method. From the vertical comparison, we can see that, as the step size decreases, the accuracy of Euler method, heun method, and Runge-Kutta method increases.

It is still necessary to make a general description of the X = x (t) and Y = Y (t) images:

At the beginning, the number of rabbits and foxes exceeded the number of foxes in the ecosystem (in other words, the number of rabbits could not meet the needs of foxes, due to lack of food, rabbits die because too many Foxes prey to reduce the number of foxes. When the number of foxes is reduced to a certain extent, rabbits continue to breed due to the decrease of natural enemies, so the number increases, at the same time, as the number of foods increases, the fox gradually begins to breed when there is enough food. As a result, the growth rate of rabbits decreases to a certain extent, the number of rabbits does not meet the needs of foxes. As a result, foxes begin to decrease ......

For Y-X images, we have already said that this is a closed curve, which illustrates the self-regulating process of the ecosystem and reveals the dynamic relationship between the number of prey and the number of prey, the close of the image also demonstrates the correctness of the model: There is a periodic oscillating Relationship Between the number displayed in the cycle: This is a long period of time, within a period of time, add together-the relationship between x (t)-T, y (t)-T is explained from another aspect.

It must be emphasized again that

The critical point in the models of the prey and the prey is (D/C, A/B) = (4000,100), and it is known that the critical point is center. The closer the Initial Value Condition x0 is (4000,100), the more the periodic solution is like the elliptic solution of the corresponding linear equations. The characteristic value of G' (4000,100) is λ = ± √ (AD) * I = ± 2/5 * √ 10i, so the period is (√ 10) π/2, or approximately 4.97.