Commit 2af4d58e authored by SOLIMAN Sylvain's avatar SOLIMAN Sylvain

invariants integrated into conservations, more or less

parent c6039c01
......@@ -66,18 +66,18 @@ check_conservations :-
forall(
(
% model: current_model? compute_ode already supposes a single model
% otherwise ofr each selected model, get all conservations and check
item([kind: conservation, id: Id]),
% otherwise for each selected model, get all conservations and check
item([kind: conservation, item: Conserv, id: Id]),
get_annotation(Id, conservation_list, C)
),
check_conservation(C)
check_conservation(C, Conserv)
).
search_conservations :-
biocham_command,
doc('calls \\command{search_conservations/1} with the default size of 4.'),
true.
search_conservations(4).
search_conservations(Integer) :-
......@@ -88,7 +88,7 @@ search_conservations(Integer) :-
size \\argument{Integer}. Such P-invariants are particular mass
conservation laws that are independent from the kinetics.'
),
true.
invariants:find_pinvar(Integer).
%%% End of public API
......@@ -100,11 +100,11 @@ solution_to_conservation(Solution, Conservation) :-
reaction_editor:simplify_solution(List, Conservation).
check_conservation(C) :-
check_conservation(C, Conserv) :-
(
check_conserv_reactions(C)
->
print_message(informational, conservation_checked(C, reactions, true))
print_message(informational, conservation_checked(Conserv, reactions, true))
;
print_message(
informational,
......@@ -112,10 +112,10 @@ check_conservation(C) :-
),
check_conserv_num(C)
->
print_message(informational, conservation_checked(C, kinetics, true))
print_message(informational, conservation_checked(Conserv, kinetics, true))
;
print_message(informational, conservation_checked(C, kinetics, false)),
print_message(warning, conservation_trusted(C))
print_message(informational, conservation_checked(Conserv, kinetics, false)),
print_message(warning, conservation_trusted(Conserv))
).
......@@ -136,7 +136,6 @@ check_conserv_reactions(C) :-
check_conserv_num(C) :-
get_weighted_ode(C, WeightedSum),
simplify(WeightedSum, Sum),
write(Sum), nl,
Sum == 0.
......
:- use_module(library(clpfd)).
:- module(
invariants,
[
find_pinvar/1
]
).
:- reexport(library(clpfd)).
:- use_module(library(lists)).
cpu_time(T) :-
statistics(walltime, [T, _]).
% avoid useless symmetrical rules
:- dynamic(invar_done/2).
% store one in, one out places/transitions
......@@ -23,12 +25,10 @@ cpu_time(T) :-
% find the complete set of minimal P-invariants
find_all_pinvar :-
biocham_command,
find_pinvar(sup).
find_pinvar :-
biocham_command,
find_pinvar(4).
......@@ -42,7 +42,6 @@ find_pinvar(ForcedMax) :-
find_invar_aux(ForcedMax, Type, OtherType, Operator) :-
reaction_graph,
cpu_time(Time1),
findall(
(Ins, Outs),
(
......@@ -70,7 +69,7 @@ find_invar_aux(ForcedMax, Type, OtherType, Operator) :-
),
% singletons
vars(Vars),
writeall(
write_list_to_solution_if(
(
Goal2 =.. [Type, P],
Goal2,
......@@ -78,10 +77,7 @@ find_invar_aux(ForcedMax, Type, OtherType, Operator) :-
assertz(base_mol([P]))
),
[P]
),
cpu_time(Time2),
Time is Time2 - Time1,
format_debug(3, "Time: ~w ms~n", [Time]).
).
arc(A, B, C) :-
......@@ -137,7 +133,6 @@ mark_singular(Type) :-
find_invar(IOList, ForcedMax, Operator) :-
retractall(invar_done(_, _)),
remove_common(IOList, UIOList),
format_debug(8, "~w no_common -> ~w~n", [IOList, UIOList]),
% add constraints X.M=0
get_constraints(UIOList, [], [], Vars, VarList, 1, MaxDomain1, Operator),
retractall(vars(_)),
......@@ -151,56 +146,41 @@ find_invar(IOList, ForcedMax, Operator) :-
),
% add a large upper bound for minimal invariants
VarList ins 0..MaxDomain,
format_debug(4, "~w~n", [MaxDomain]),
format_debug(4, "~w~n~w~n", [Vars, VarList]),
build_sum(VarList, Sum),
format_debug(8, "Sum: ~w~n", [Sum]),
!,
% avoid the trivial solution
Sum #> 0,
format_debug(9, "Sum #> 0~n", [Sum]),
nb_setval(base, []),
format_debug(9, "looking for paths~n", []),
find_paths,
format_debug(9, "found paths~n", []),
% zero the alternative of each path
zero_paths(Vars, VarList),
format_debug(9, "zeroed paths~n", []),
% reorder vars by neighborhood
% FIXME compare with reverse(VarList, ListVar),
reorder(Vars, VarList, ListVar),
format_debug(9, "reordered vars~n", []),
(
repeat,
(
get_base_constr(VarList),
cpu_time(Time1),
% does NOT guarantee a base e.g.
% add_rule(2*Y + Z => 3*X). for labeling not ff
% my_labeling(VarList)
label(ListVar)
->
cpu_time(Time2),
nb_getval(base, B),
(
MaxDomain > 1
->
make_base(B, VarList, BB)
;
BB = [VarList | B],
format_debug(6, "New minimal invariant: ~w~n", [VarList])
BB = [VarList | B]
),
Time is Time2 - Time1,
nb_setval(base, BB),
length(BB, LBB),
format_debug(2, "new base length: ~w~n", [LBB]),
format_debug(2, "labeling time: ~w~n", [Time]),
fail
;
retractall(base_mol(_)),
expand_base(Vars),
nb_setval(nb, 0),
writeall(
write_list_to_solution_if(
(
base_mol(Mol),
nb_getval(nb, N),
......@@ -225,24 +205,19 @@ expand_base(Vars) :-
nb_getval(base, B),
member(VL, B),
to_mols(VL, Vars, Mol),
format_debug(4, "Mol ~w~n", [Mol]),
assertz(base_mol(Mol)),
fail.
expand_base(_) :-
normalized_path(E, L),
format_debug(5, "Path ~w ~w~n", [E, L]),
findall(
[E | MMol],
(
base_mol(Mol),
format_debug(5, "applying ~w to ~w~n", [L, Mol]),
subtract(Mol, L, MMol),
format_debug(5, "-> ~w~n", [MMol])
subtract(Mol, L, MMol)
),
NewBases
),
format_debug(5, "Newbases ~w~n", [NewBases]),
assert_to_base(NewBases),
fail.
......@@ -293,7 +268,6 @@ find_paths :-
find_paths(Length, _Maxlength) :-
format_debug(8, "looking for paths of length ~w~n", [Length]),
singular(E),
arc(Before, E, 1),
arc(E, After, 1),
......@@ -302,7 +276,6 @@ find_paths(Length, _Maxlength) :-
->
retract(singular(E)),
assertz(path(E, Path)),
format_debug(3, "~w -> ~w~n", [E, Path]),
fail
).
......@@ -380,9 +353,6 @@ get_constraints([(KL, KR) | L], Vars, VarList, NVars, NVarList, MaxDomain, NMaxD
Constr =.. [Operator, SL, SR],
Constr,
% add constraint X.M = O
format_debug(3, "~w~n", [Constr]),
format_debug(3, "Max:~w = ~w * max(~w, ~w)~n~w~n",
[MaxDomain2, MaxDomain, Max1, Max2, Constr]),
get_constraints(L, Vars2, VarList2, NVars, NVarList, MaxDomain2,
NMaxDomain, Operator).
......@@ -461,9 +431,7 @@ to_mols([N | T], [V | Vars], [N * V | Mols]) :-
get_base_constr(VarList) :-
nb_getval(base, B),
format_debug(3, "Base: ~w~n", [B]),
get_base_constr_rec(B, VarList),
format_debug(3, "got constr...~n", []).
get_base_constr_rec(B, VarList).
get_base_constr_rec([], _).
......@@ -496,7 +464,6 @@ make_base(T, H, [H | T]) :-
% if new vector has max value 1, cannot subsume an earlier vect
% FIXME if MaxDomain=1 this check is useless...
max_list(H, 1),
format_debug(5, "skipping subsumption test for ~w~n", [H]),
!.
make_base([], V, [V]).
......@@ -505,7 +472,6 @@ make_base([H | T], V, B) :-
(
bigger_support(H, V)
->
format_debug(3, "removing a vector from the base~n", []),
make_base(T, V, B)
;
B = [H | B1],
......@@ -573,7 +539,6 @@ remove_common([(SL, M) | KL], KR, UKL, UKR) :-
% reorders VarList1 into VarList2 s.t. most neighbors (computed using Values)
% are taken first
reorder(Vars, L1, L2) :-
format_debug(8, "Original order: ~w~n", [Vars]),
% initialize neighbors to 0 for everyone
make_keylist(Vars, KV),
reorder_aux(KV, L1, L2).
......@@ -588,7 +553,6 @@ reorder_aux(KL, VL, [V | L]) :-
keysort(KL, SKL),
% sort is ascending, so biggest is last...
last(SKL, _-H),
format_debug(9, "~w ", [H]),
nth1(N, KL, _-H),
% get corresponding variable
nth1(N, VL, V),
......@@ -633,52 +597,18 @@ update_neighbors([N-A | L], B, [M-A | LL]) :-
),
update_neighbors(L, B, LL).
%% format_debug(+Level, +String, +Args)
%
% print debug messages if debug above Level
format_debug(Level, String, Args):-
catch(
nb_getval(debug, N),
error(existence_error(variable, debug), _),
N = 0
),
(
(N >= Level)
->
write_list_to_solution_if(Goal, List) :-
forall(
Goal,
(
(
(Level > 0)
->
write('** Debug ** ')
;
true
),
format(String, Args),
flush_output
reaction_editor:list_to_solution(List, Solution),
write_term(Solution, [quoted(false)]),
nl
)
;
true
).
%% writeall(?Goal, ?Term)
%% writeall(+Stream, ?Goal, ?Term)
%
% calls Goal and writes to Stream all possible Terms as instantiated by Goal
writeall(Goal, Term) :-
writeall(user_output, Goal, Term).
writeall(Stream, Goal, Term) :-
call(Goal),
write_term(Stream, Term,[quoted(false)]),
nl(Stream),
fail.
writeall(_, _, _).
%% build_sum(+List, ?Sum)
%
% builds the sum of a list
......
Markdown is supported
0% or
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment