Stability of fixed points in the Blanchard, Devaney, Hall model of trees, moose, and wolves.

f[x_, y_, z_] = r[1]*x*(1 - x) - a[1]*x*y

g[x_, y_, z_] = r[2]*y*(1 - y) + a[2]*x*y - a[3]*y*z

h[x_, y_, z_] = r[3]*z*(1 - z) + a[4]*y*z

sol = 
      Simplify[Solve[{f[x, y, z] == 0, g[x, y, z] == 0, 
            h[x, y, z] == 0}, {x, y, z}]]

sol[[5]]

Jac = {{D[f[x, y, z], x], D[f[x, y, z], y], 
          D[f[x, y, z], z]}, {D[g[x, y, z], x], D[g[x, y, z], y], 
          D[g[x, y, z], z]}, {D[h[x, y, z], x], D[h[x, y, z], y], 
          D[h[x, y, z], z]}}

MatrixForm[Jac]

At this point, the formulas get completely out of hand unless we pick some values for the various parameters.

par = {r[1] -> 1, r[2] -> 2, r[3] -> 2, 
        a[1] -> 0.5, a[2] -> 1, a[3] -> 1, 
        a[4] -> 1}

JJ = Jac /. par

Now evaluate the Jacobian matrix at the fixed point.

JJstar = JJ /. {sol[[5]] /. par}

Finally compute the eigenvalues of the Jacobian.

N[Eigenvalues[JJstar]]

Thus, this population is stable.