**Background**

To use permanence, Lotka-Volterra dynamics have to be assumed, because it is only in this case that a sufficient condition for permanence is known:

Such a dynamical system is permanent if two conditions hold (Hofbauer and Sigmund 1988, *The theory of evolution and dynamical systems* p. 98). (1) It is dissipative; this is true for Lotka-Volterra systems where all the basal species are self-limiting and heterotrophs cannot survive without their prey. (2) That:

for some is an average Lyapunov function.

In dissipative Lotka-Volterra systems, the test for permanence reduces to a linear programming problem (Jansen 1987, *J. Math. Biol*.; Law & Morton 1996, *Ecology*).

An alternative way of describing the system, pursued in an unpublished paper by Richard Law and R. Daniel Morton, is in terms of the convex hull and the set where the densities of all species are non-increasing. The convex hull of the boundary equilibria is the smallest convex set of which every boundary equilibrium is a member. If this set is disjoint from , then the system is permanent. So to test for permanence, one can find a hyperplane passing through the interior equilibrium that has a greater value than that of the hyperplane passing through every boundary equilibrium. One can also measure the ‘strength’ of permanence by , which is the shortest distance from to , which is always the shortest distance from to the interior equilibrium (Figure 1).

**Example**

You can download the complete example for Octave from my website: mortongen.m.

Consider the four-species system shown in Figure 2 with the parameter values below, adapted from an unpublished paper by Richard Law and R. Daniel Morton, with:

and

It has an interior equilibrium point:

The first step is to find all subsystems , and then identify which of these have a positive boundary equilibrium .

If we had all subsystems in `subStore`

, then equilibria would be evaluated by:

for subCnt = 2:noSubs; % Acquire list of species absent and present present = subStore(subCnt,:); present(find(present == 0)) = []; absent = 1:nospp; absent(present) = []; % Evaluate boundary equilibrium ANow = A; ANow(absent,:) = []; ANow(:,absent) = []; rNow = r; rNow(absent) = []; if rank(ANow) == size(ANow,1); xSub = ANow-rNow; % etc.

For our example system, the positive boundary equilibria are:

The transversal eigenvalues at the boundary equilibria

are found with code like:

rowofG = A*(xStore(subCnt,:)')+r;

and are:

The linear programming problem is then to minimise subject to:

where each positive equilibrium point on the boundary gives rise to a constraint, and is part of the average Lyapunov function. So if one can find a solution with , then one can be certain that is an average Lyapunov function.

An introduction on how to code linear programming problems in Octave can be found in a previous post. The vector for the Lyapunov function, as found using Jansen’s linear program method, is

Hofbauer & Sigmund (1988, p. 176-177) states that the convex hull can be separated from by a hyperplane with a “separating functional” , where is the same used in the Lyapunov function. This separating functional is the vector normal to the hyperplane.

The equation for the hyperplane is

If one wants to go straight to finding , then substitute in the boundary steady state for each positive subsystem, reducing the value of until the hyperplane is above or equal to every boundary steady state.

So for subsystem ,

which gives

Then subsystem is next, with

Substituting into the hyperplane equation gives

which implies that this boundary steady state is above the hyperplane, therefore it is used to find a revised value of ,

This process is continued, checking that all other boundary steady

states are below the plane, and revising if they are not, which gives the final value of .

Now the distance from the interior equilibrium to the plane can be

found. The interior equilibrium is

and so the perpendicular distance from the hyperplane to is

Alternatively can be calculated for each boundary steady state, and the smallest taken as the final value. If `xStore`

stores the values of all the boundary equilibria, then:

d = (n'*xs-xStore(posSubs,:)*n)/norm(n); disp('d = ') disp(d) mind = min(d); disp('smallest distance = ') disp(mind)

—

You can download the complete example for Octave from my website: mortongen.m.

Hey, I just stopped in to visit your website and thought I’d say I enjoyed myself.

Thank you Jailbreak, I rather enjoy measuring the distance from the convex hull to the interior equilibrium myself.