CHW Example 6.6

Loads the chemJac command file, which contains all commands used below

<< "chemJac.m"

Inputting the stoichiometric matrix S

S = {{-1, 1,  -1,  1,  0,  0,  0,  0}, {-1, 1,  0,  
          0,  1,  -1,  0,  0},  {0,  0,  -1,  1,  -1,  
          1,  -1,  1}, {0,  0,  0,  0,  -1,  1,  
          1,  -1}, {0,  0,  0,  0,  0,  0,  -1,  
          1}, {1,  -1,  0,  0,  0,  0,  0,  0}, {0,  0,  
          1,  -1,  0,  0,  0,  0}};  S // MatrixForm

Making the non-mass-action flux vector 'monomials' and the list of species concentrations 'vars' automatically from S

monomials = makeMonomialNMA[S]

vars = svars[S]

'S.monomials' is the right hand side of the ODE

S . monomials

The Jacobian of the RHS of the ODE, and its Craciun-Feinberg determinant

j = jac[S . monomials,  vars];
j // MatrixForm

det = cfDet[j];

coeffs[det]

(k[6])'[a[2]] (k[8])'[a[4]] Derivative[0, 1][k[3]][a[1], a[3]] Derivative[1, 0][k[1]][a[1], a[2]];

Collect[det, (k[8])'[a[4]] Derivative[0, 1][k[3]][a[1], a[3]] Derivative[1, 0][k[1]][a[1], a[2]]] ;

(-1+(k[6])'[a[2]]) (k[8])'[a[4]] Derivative[0, 1][k[3]][a[1], a[3]] Derivative[1, 0][k[1]][a[1], a[2]];(*term resulting from the previous command ...  we see that as long as (-1+(k[6])'[a[2]]) is not positive, the determinant will be strictly negative on the open orthant of concentration space*)

Collect[det,Derivative[0, 1][k[3]][a[1], a[3]] Derivative[1, 0][k[1]][a[1], a[2]]] ;

Derivative[0, 1][k[3]][a[1], a[3]] (-1-(k[8])'[a[4]]+(k[6])'[a[2]] (k[8])'[a[4]]-Derivative[0, 1][k[5]][a[3], a[4]]-Derivative[0, 1][k[7]][a[3], a[5]]-Derivative[0, 1][k[5]][a[3], a[4]] Derivative[0, 1][k[7]][a[3], a[5]]) Derivative[1, 0][k[1]][a[1], a[2]];

Simplify[-1-(k[8])'[a[4]]+(k[6])'[a[2]] (k[8])'[a[4]]-Derivative[0, 1][k[5]][a[3], a[4]]-Derivative[0, 1][k[7]][a[3], a[5]]-Derivative[0, 1][k[5]][a[3], a[4]] Derivative[0, 1][k[7]][a[3], a[5]]]

Core Determinant of the Jacobian of the RHS of the ODE

core = coreDet[j, S]

coeffs[core]