## Description

We can also use the RK4 method to numerically solve systems of differential equations. Consider the predator-prey system

(1)As a nonlinear system, it is difficult, if not impossible, to find an analytical solution. However, we can use `desolve_system_rk4` to approximate the solution curves. Since the command only outputs a list of points, we'll have to use a bit of Python to construct the plot ourselves. The arguments of the command are pretty much the same as those for `desolve_system`, except the command only takes the right-hand side of each equation and there are additional arguments `end_points`, which specifies the range over which we are approximating, and `step`, which specifies the step size for the RK4 method. A lower step size increases precision, but also computation time:good options are `0.1` and `0.01`, depending on your approximation range.

In the first cell, we'll plot both our solution curves with respect to $t$:

## Sage Cell

#### Code

```
var('x y t')
dx = 2*x - x*y
dy = -5*y + x*y
pts = desolve_system_rk4(des=[dx, dy], vars=[x, y], ivar=t, step = 0.01, end_points=[0, 10], ics=[0, 1, 1])
ptsx = [[i, j] for i, j, k in pts]
ptsy = [[i, k] for i, j, k in pts]
p = line(ptsx, color='blue') + line(ptsy, color='red')
p.show(axes_labels=['$t$', '$x, y$'])
```

In this cell, we plot the solution as a parametric curve in the $xy$ plane:

#### Code

```
var('x y t')
dx = 2*x - x*y
dy = -5*y + x*y
pts = desolve_system_rk4(des=[dx, dy], vars=[x, y], ivar=t, step = 0.01, end_points=[0, 10], ics=[0, 1, 1])
ptspar = [[j, k] for i, j, k in pts]
p = line(ptspar, color='blue')
p.show(axes_labels=['$x$', '$y$'], xmin=0, ymin=0)
```

In this cell, we do the same as in the previous cell, but superimpose the curve over the vector field of the system:

#### Code

```
var('x y t')
dx = 2*x - x*y
dy = -5*y + x*y
mag = sqrt(dx^2 + dy^2)
vf = vector([dx, dy])
pts = desolve_system_rk4(des=[dx, dy], vars=[x, y], ivar=t, step = 0.01, end_points=[0, 10], ics=[0, 1, 1])
ptspar = [[j, k] for i, j, k in pts]
p = line(ptspar, color='blue') + plot_vector_field(vf/mag, (x, 0, 15), (y, 0, 10))
p.show(axes_labels=['$x$', '$y$'], xmin=0, ymin=0)
```

## Options

none

## Tags

Primary Tags:

Secondary Tags:

A list of possible tags can be found at The WeBWorK Open Problem Library. For linear algebra tags see the Curated Courses Project.

## Related Cells

Any related cells go here. Provide a link to the page containing the information about the cell.

## Attribute

Permalink:

Author: Thomas Judson

Date: 14 May 2020 05:06

Submitted by: Zane Corbiere