How should I modify the rate rule using Symbiology so that it works?

I am attempting to model a system for the Guerbert reaction. The reaction begins with Ethanol, which couples with itself to form Butanol, which couples with Ethanol to form Hexanol....and so on.

Here's the code:

 m1 =  sbiomodel('Guerbert Reaction');

c1 = addcompartment(m1, 'Plug-flow reactor');
% compartment added

s1 = addspecies(m1, 'EtOH');
s2 = addspecies(m1, 'BtOH');
s3 = addspecies(m1, 'HxOH');
% create all species

set(s1, 'InitialAmount', 5)
set(s2, 'InitialAmount', 0)
set(s3, 'InitialAmount', 0)
%set initial amounts for all species

r1 = addreaction(m1, 'EtOH + EtOH -> BtOH');
r2 = addreaction(m1, 'EtOH + BtOH -> HxOH');
%reactions added, parent object is m1 model

p1 = addparameter(m1, 'k2_2', 'Value', 0.1);
p2 = addparameter(m1, 'k2_4', 'Value', 0.4);
%parameters added to model object. This ensures that all objects within the model can use it

k1 = addkineticlaw(r1, 'MassAction');
k2 = addkineticlaw(r2, 'MassAction');
% kinetic laws added, of Mass Action type. Parent is the reaction object

p1 = addparameter(k1, 'k2_2', 'Value', 0.1);
p2 = addparameter(k2, 'k2_4', 'Value', 0.4);
%parameters added for use within the kinetic law

setparameter(k1, 'Forward Rate Parameter', 'k2_2')
setparameter(k2, 'Forward Rate Parameter', 'k2_4')
% defining the kinetic law by specifying which parameters will be substituted in the expression

ru1 = addrule(m1, 'EtOH = (-k2_2 * EtOH)/(1+EtOH) - (k2_4 * BtOH)/(1+EtOH)')
ru2 = addrule(m1, 'BtOH = (k2_2 * EtOH)/(1+EtOH) - (k2_4*BtOH)/(1+EtOH)')
ru3 = addrule(m1, 'HxOH = (k2_4 * BtOH)/(1*EtOH)')
% rate rules

set(ru1, 'RuleType', 'Rate')
set(ru2, 'RuleType', 'Rate')
set(ru3, 'RuleType', 'Rate')
% rule type is "rate" because the rules are differential}

However, when I verify the model, I get this error message:

>> verify(m1)
Error using SimBiology.Model/verify
--> Error reported from ODE Compilation:
Invalid rule variable 'EtOH' in rate rule 'EtOH = (-k2_2 * EtOH)/(1+EtOH) - (k2_4 * BtOH)/(1+EtOH)'. This object is being used in a reaction or in another rule.

What does this mean? How may I correct this?

Answer 1

[Jeremy Huard]

This is an old question but an answer could still benefit others who bumped into the same issue.

The key piece of the error message is: This object is being used in a reaction or in another rule

Indeed, the ODE driving the dynamics of a species can only be defined by EITHER a set of reactions OR a rate rule. But not both.

For the code above to work, 3 options can be considered:

1. remove the reactions from the model and only keep the rate rules
2. remove the rate rule and only keep the reactions. There one can use a built-in kinetic law or type in a custom reaction rate.
3. one could keep the reactions for visualization purposes and the rate rules. but for this to work, the species property BoundaryCondition must be set to true.

Side note:

There is no need to add the parameters to the kinetic laws (lines of code with addparameter(kx,...) ). The forward rate parameters can be the model-scoped parameters.

The code above creates two independent k2_2 parameters: one scoped to the model (accessible by all reactions), one scoped to the reaction (only accessible by this reaction). Changing one does not change the other.

Here, both lines of code with addparameter(kx,...) can be removed.