"""
Ported by SCL, 27 Oct 2012.
...from combination of predprey.m and predprey_rh.m
The original top-of-file comments appear in function docstrings.
"""
import numpy as np
def predprey_rh(x, t, rnom, predprey_K, predprey_xd, predprey_ud):
"""
predprey_rh.m - feedback function for predator prey system
RMM, 13 Sep 06
This function computes a simple linear control law, with saturation, for
the predator prey system
"""
# Compute the feedback control law
rh = (-np.dot(predprey_K, (x - predprey_xd)) + predprey_ud)[0]
# Saturate the control input between 0 and 4*rnom
if rh < 0:
rh = 0
if rh > 4*rnom:
rh = 4*rnom
return rh
def predprey(x, t,
r = 1.6, d = 0.56, b = 0.6, k = 125, a = 3.2, c = 50,
control_args=None):
"""
predprey.m - updated predator prey model for CDS 101/110
RMM, Oct 07
This file contains a simple model for a predator prey system. The
model is developed following the treatment in the book by Murray
(the other one).
The states for these equations are
x(1) # of rabbits (prey)
x(2) # of foxes (predator)
The equations allow for a quadratic term in both the predator and
prey species that cause the population to fall off if you get too
many of a given species. These terms are turned off by default.
"""
# Extract the states
x1 = x[0]
x2 = x[1]
# Check to see if a feedback law has been defined (when r < 0)
if r < 0 and (control_args is not None):
# Assume that the user has defined a function that we can call
# Add this value to the nominal value of the parameter (renegated)
r = -r + predprey_rh(x, t, -r, *control_args)
dx1 = r * x1 * (1 - x1/k) - a * x2 * x1/(c + x1)
dx2 = b *a *x2 * x1 / (c + x1) - d * x2
return np.array([dx1, dx2])