canonical representation of solutions in reactions

rate_independence.pl

@@ -54,8 +54,6 @@ test_rate_independence_inputs_sinks :-
   doc('Test graphical sufficient conditions for rate independence of the current model for the computation of graph output species from graph input species (assuming well-formed kinetics). The input species are the species that are not a reaction product or a strict catalyst (i.e. with same stoichiometry). The output species are the species that are not reactant.'),
   input_species(Inputs),
   sink_species(Outputs),
-  write('inputs: '), writeln(Inputs),
-  write('outputs: '), writeln(Outputs),
   test_rate_independence(Inputs, Outputs).
 
 
@@ -65,6 +63,8 @@ rate_independence_reduction(Inputs, Outputs):-
   type(Inputs, {object}),
   type(Outputs, {object}),
   doc('Reduces the current reaction model while preserving rate independence for computing \\argument{Outputs} species from \\argument{Inputs} species.'),
+  write('inputs: '), writeln(Inputs),
+  write('outputs: '), writeln(Outputs),
   simplify(Inputs, Outputs),
   !,
   writeln('Reduced'),
@@ -167,7 +167,7 @@ dfs_no_loop(GraphId, Vertex) :-
 
 
 
-% Elimination of species with no outgoing fork and trivial loops
+% Elimination of trivial reactions and isolated species with no outgoing fork
 simplify(Inputs, Outputs):-
   \+ trivial_reaction,
   \+ same_reaction,

rate_independence_test.bc

+
+prolog('writeln("%%%%%%%% WRONG ANSWER %%%%%%%%%")').
+
+
+clear_model.
+d=>a.
+a=>b.
+b=>2*d.
+a=>c.
+list_model.
+test_rate_independence_inputs_sinks.
+prolog('writeln("ANSWER NO AAAARGH NON CONFLUENCE OF RESULT / REACTION ORDERING AS SHOWN BELOW:")'). 
+
+
+clear_model.
+a=>c.
+a=>b.
+b=>2*d.
+d=>a.
+list_model.
+test_rate_independence_inputs_sinks.
+prolog('writeln("ANSWER NO ")'). prolog('writeln("============")').
+
+
+
 prolog('writeln("%%%%%%%% DIVERGING SYLVAIN %%%%%%%%%")').
 
+
 clear_model.
 a=>b.
 b=> a+c.
 list_model.
 test_rate_independence_inputs_sinks.
-prolog('writeln("ANSWER DIV ")'). prolog('writeln("============")').
+prolog('writeln("ANSWER YES DIV ")'). prolog('writeln("============")').
 
 
 clear_model.
@@ -14,34 +40,44 @@ b=> a+c.
 d=> a.
 list_model.
 test_rate_independence_inputs_sinks.
-prolog('writeln("ANSWER DIV ")'). prolog('writeln("============")').
+prolog('writeln("ANSWER YES DIV ")'). prolog('writeln("============")').
 
 
 
 clear_model.
+a=>c.
 a=>b.
 b=>2*a.
-a=>c.
 list_model.
 test_rate_independence_inputs_sinks.
-prolog('writeln("ANSWER DIV ")'). prolog('writeln("============")').
+prolog('writeln("ANSWER YES DIV ")'). prolog('writeln("============")').
 
 
 clear_model.
+a=>c.
 a=>b.
-b=>c.
-c=>2*a.
-b=>d.
+b=>2*d.
+d=>a.
 list_model.
 test_rate_independence_inputs_sinks.
-prolog('writeln("ANSWER DIV ")'). prolog('writeln("============")').
+prolog('writeln("ANSWER YES DIV ")'). prolog('writeln("============")').
 
 
-prolog('writeln("%%%%%%%% WRONG ANSWER %%%%%%%%%")').
 
+prolog('writeln("============")').
 prolog('writeln("%%%%%%%% UNDECIDED %%%%%%%%%")').
 
 
+clear_model.
+a => x+c.
+b => y+c.
+x+y => z.
+c+z => r.
+list_model.
+test_rate_independence_inputs_sinks.
+prolog('writeln("ANSWER YES: max = sum-min")'). prolog('writeln("============")').
+
+
 clear_model.
 a+b=>c+d.
 c+d=>a+b.
@@ -52,17 +88,26 @@ test_rate_independence_inputs_sinks.
 prolog('writeln("ANSWER YES: composite loop leakage with and-fork")'). prolog('writeln("============")').
 
 
+prolog('writeln("%%%%%%% CORRECT ANSWER %%%%%%%%")').
+
+
+clear_model.
+2*a => a+b.
+a+b => 2*a.
+a+b => 2*b.
+list_model.
+test_rate_independence_inputs_sinks.
+prolog('writeln("ANSWER NO: subtle multiple steady states ")'). prolog('writeln("============")').
+
+
 clear_model.
 MA(k1) for a+b => b.
 MA(k2) for b+c=>d.
 list_model.
 test_rate_independence_inputs_sinks. % no
-prolog('writeln("ANSWER YES: irrelevant input")'). prolog('writeln("============")').
-
+prolog('writeln("ANSWER NO: for on input a")'). prolog('writeln("============")').
 
 
-prolog('writeln("%%%%%%% CORRECT ANSWER %%%%%%%%")').
-
 
 clear_model.
 MA(k1) for a+b => b.
@@ -92,16 +137,6 @@ test_rate_independence_inputs_sinks.
 prolog('writeln("ANSWER YES: composite double loop leakage with composite or-fork")'). prolog('writeln("============")').
 
 
-clear_model.
-a => x+c.
-b => y+c.
-x+y => z.
-c+z => r.
-list_model.
-test_rate_independence_inputs_sinks.
-prolog('writeln("ANSWER YES: max = sum-min")'). prolog('writeln("============")').
-
-
 clear_model.
 a => b.
 b => c.

reaction_editor.pl

@@ -689,10 +689,21 @@ make_reaction(
   ).
 
 
+:- devdoc(' The left and right hand sides of reactions are now sorted (canonical form) and contain the catalysts').
 
 make_reaction(
-  Name, Kinetics, Left, Inhibitors, Catalyst, Right, Reversible, ReactionNamed
-) :-
+  Name, Kinetics, UnsortedLeft, UnsortedInhibitors, UnsortedCatalyst, UnsortedRight, Reversible, ReactionNamed
+	     ) :-
+  append(UnsortedLeft, UnsortedCatalyst, UnsortedLC),
+  simplify_solution(UnsortedLC,SimplifiedLC),
+  
+  append(UnsortedRight, UnsortedCatalyst, UnsortedRC),
+  simplify_solution(UnsortedRC,SimplifiedRC),
+  
+  sort(SimplifiedLC, Left),
+  sort(SimplifiedRC, Right),
+  sort(UnsortedInhibitors, Inhibitors),
+ 
     (
         Name=''
      ->
@@ -703,7 +714,6 @@ make_reaction(
   Reaction = (Kinetics for Body),
   list_to_solution(Left, LeftSolution),
   list_enumeration(Inhibitors, InhibitorsEnumeration),
-  list_to_solution(Catalyst, CatalystSolution),
   list_to_solution(Right, RightSolution),
   (
     InhibitorsEnumeration = '_'
@@ -715,25 +725,9 @@ make_reaction(
   (
     Reversible = true
   ->
-    (
-      CatalystSolution = '_'
-    ->
-      Body = (Reactants <=> RightSolution)
-    ;
-      simplify(Reactants + CatalystSolution, Reactants_full),
-      simplify(RightSolution + CatalystSolution, Products_full),
-      Body = (Reactants_full <=> Products_full)
-    )
+   Body = (Reactants <=> RightSolution)
   ;
-    (
-      CatalystSolution = '_'
-    ->
-      Body = (Reactants => RightSolution)
-    ;
-      simplify(Reactants + CatalystSolution, Reactants_full),
-      simplify(RightSolution + CatalystSolution, Products_full),
-      Body = (Reactants_full => Products_full)
-    )
+   Body = (Reactants => RightSolution)
   ).
 
 

reaction_editor.plt

@@ -14,7 +14,7 @@ test(
 
 test(
   'catalyst',
-  [Reactions == ['MA'(1) for 2*b+a+c => b+a+c, 'MA'(1) for b+a+c => 2*b+a+c]]
+  [Reactions == ['MA'(1) for a+c+2*b => a+b+c, 'MA'(1) for a+b+c => a+c+2*b]]
 ) :-
   clear_model,
   add_reaction(a + b + c <=[ b ]=> a + c),
@@ -31,7 +31,7 @@ test(
 
 test(
 'michaelis menten reduction',
-  [Reactions == ['E'*'S'*'C' for 'S'+'E'+'C'=> 'P'+'E'+'C']]
+  [Reactions == ['E'*'S'*'C' for 'C'+'E'+'S'=> 'C'+'E'+'P']]
 ) :-
   command(load('library:examples/michaelis-menten/mm.bc')),
   command(delete_reaction('C'=>'E'+'S')),
@@ -40,7 +40,7 @@ test(
 
 test(
 'michaelis menten reduction C',
-  [Reactions == ['E'*'S'*'C' for 'S'+'E'+'C'=> 'P'+'E'+'C']]
+  [Reactions == ['E'*'S'*'C' for 'C'+'E'+'S'=> 'C'+'E'+'P']]
 ) :-
 
   command(load('library:examples/michaelis-menten/mm.bc')),

reaction_graphs.plt

@@ -24,7 +24,7 @@ test(
 
 test(
   'import_reactions_from_graph',
-  [true(Reactions == ['MA'(1) for a+b=>c+b, 'MA'(2) for a + 3 * c => 2 * d ])]
+  [true(Reactions == ['MA'(1) for a+b=>b+c, 'MA'(2) for a + 3 * c => 2 * d ])]
 ) :-
   clear_model,
   new_graph,