Commit 7cea27c7 authored by Mathieu Hemery's avatar Mathieu Hemery
Browse files

Reflow quite a lot of list processing that was done in imperative way to natural prolog

parent c4fec086
......@@ -177,67 +177,48 @@ hybrid_static_simulation(ODEFilename, StochFilename, OutFilename, Volume, Time)
build_constraints_by_rules(Reactants_by_rules, SortedReactants, 1, [], Hybrid_reactants_list),
build_constraints_list(Hybrid_reactants_list, Constraints_list),
findall(Reactant_stoch,
findall(Reactant_pre,
(
member(Reaction,Reactions2),
reaction(Reaction, [kinetics: _,reactants: Reactants,inhibitors: _, products: _]),
member(Reactant,Reactants),
Reactant = _*Reactant_pre,
species_to_stoch(Reactant_pre, Reactant_stoch)
Reactant = _*Reactant_pre
),
Reactant_list),
findall(Sp_stoch,
(
member(Sp, Reactant_list),
species_to_stoch(Sp, Sp_stoch)
),
Reactant_stoch_list),
findall(Product_stoch,
findall(Product_pre,
(
member(Reaction,Reactions2),
reaction(Reaction, [kinetics: _,reactants: _,inhibitors: _, products: Products]),
member(Product,Products),
Product = _*Product_pre,
species_to_stoch(Product_pre, Product_stoch)
Product = _*Product_pre
),
Product_stoch_list),
append(Reactant_stoch_list, Product_stoch_list, StochSpecies),
sort(StochSpecies,SortedStochSpecies), % To eliminate duplication.
% build StochReactants(with '_total') from SortedStochSpecies.
length(SortedStochSpecies, Species_number),
build_reactants(SortedStochSpecies, SortedODESpecies, Species_number, 1, [], StochReactants),
build_reactants_old(SortedStochSpecies, SortedODESpecies, Species_number, 1, [], StochReactantsold),
writeln(SortedStochSpecies),
writeln(SortedODESpecies),
writeln(StochReactants),
writeln(StochReactantsold),
%========================================%
%==== Build Data Lists for Reactions ====%
%========================================%
prepare_change_list(Reactions2,ODESpecies,List),
findall([Species_change, Alpha],
Product_list),
findall(Sp_stoch,
(
member(Element,List),
% To get [RSpecies, PSpecies] from [Reactants, Products]
Element = [Reactants,Products,Alpha],
findall([Species1, Minus_number],
(
member(R, Reactants),
R = Minus_number * Rname,
species_to_stoch(Rname, Species1)
),
RSpecies),
findall([Species2, Add_number],
(
member(P, Products),
P = Add_number * Pname,
species_to_stoch(Pname, Species2)
),
PSpecies),
%== Build list (Species_change) to record species change ==%
append_reactants(RSpecies, 1, [], Temp_list, SortedStochSpecies),
append_products(PSpecies, 1, Temp_list, Changes, SortedStochSpecies),
build_change_list(Changes, 1, Species_number, [], Species_change)
member(Sp, Product_list),
species_to_stoch(Sp, Sp_stoch)
),
Reactions_data),
Product_stoch_list),
append(Reactant_list, Product_list, Species_list),
sort(Species_list, Sorted_Species), % To eliminate duplication.
append(Reactant_stoch_list, Product_stoch_list, StochSpecies),
sort(StochSpecies, SortedStochSpecies), % To eliminate duplication.
% build StochReactants(with '_total' or '_stoch') from Sorted_Species.
build_reactants(Sorted_Species, SortedODESpecies, StochReactants),
%== Build Events and Calculate Alphas==%
print_event(Reactions_data, SortedStochSpecies, StochReactants,Hybrid_reactants_list, Constraints_list, Stream),
%===================================================================%
%==== Build Data Lists for Reactions and print events in Stream ====%
%===================================================================%
prepare_change_list(Reactions2,ODESpecies,List),
build_reaction_data(List, Sorted_Species, Reaction_data_list),
print_event(Reaction_data_list, SortedStochSpecies, StochReactants,Hybrid_reactants_list, Constraints_list, Stream),
format(Stream, "~n% Plot part~n~n", []),
format(Stream,"option(show: {}). % You can choose here which species you want to plot : option(show: {a_total,b_total,c_total,d_total}) to plot species a,b,c and d.~n",[]),
......@@ -344,10 +325,9 @@ build_constraints_list([List| InTail], Counter, [Constraint| OutTail]) :-
build_constraints_list(InTail, Counterp, OutTail).
% build_reactants(StochSpecies, ODESpecies, Species_number, Counter, Init, StochReactants)
build_reactants(StochSpecies, ODESpecies, _SpN, _Cou, _Ini, StochReactants) :-
build_reactants(StochSpecies, ODESpecies, StochReactants).
% build_reactants(Species, ODESpecies, Reactants)
%
% Construct the list of reactant species in an adequate format
build_reactants([], _ODESpecies, []).
......@@ -362,28 +342,10 @@ build_reactants([Species|SpTail], ODESpecies, [Reactant|ReacTail]) :-
build_reactants(SpTail, ODESpecies, ReacTail).
build_reactants_old(StochSpecies, ODESpecies, Species_number, Counter, Init, StochReactants) :-
nth1(Counter, StochSpecies, This_species),
(
member(This_species, ODESpecies)
->
species_to_total(This_species, This_reactant)
;
species_to_stoch(This_species, This_reactant)
),
append(Init, [This_reactant], Next),
(
Counter < Species_number
->
Next_counter is Counter + 1,
build_reactants_old(StochSpecies, ODESpecies, Species_number, Next_counter, Next, StochReactants)
;
StochReactants = Next
).
% write_hybrid_ode(+Stream, +Hybrid_list)
write_hybrid_ode(Stream,Hybrid_original_list) :-
list_hybrid_ode(Stream,Hybrid_original_list),
write_hybrid_ode(Stream, Hybrid_original_list) :-
list_hybrid_ode(Stream, Hybrid_original_list),
with_output_to(Stream, models:list_model_initial_state),
with_output_to(Stream, models:list_model_parameters),
with_output_to(Stream, models:list_model_events),
......@@ -526,50 +488,43 @@ convert_into_hybrid2(Molecule,Hybrid_original_list,HybridMolecule) :-
atom_concat(Molecule5,'))',HybridMolecule),
!.
build_change_list(Changes, Counter, Length, Init, Result) :-
% build_stoechiometry(+Reactant_list, +Product_list, +Species_list, -Stoechiometry)
%
% Construct the stoechiometry of a reaction
% stoech are typed as atom to explicitely display a '+/-' in the output file
build_stoechiometry(_Reactant_list, _Product_list, [], []) :- !.
Counter =< Length
->
(
build_stoechiometry(Reactant_list, Product_list, [Sp|TailSp], [Stoech|TailStoech]) :-
findall(CR, member(CR*Sp,Reactant_list), CR_list),
findall(CP, member(CP*Sp,Product_list), CP_list),
sum_list(CR_list, CRT),
sum_list(CP_list, CPT),
Stoech_number is CPT-CRT,
(
\+ member([Counter, _], Changes)
Stoech_number < 0
->
(
atom_number(Zero,0),
atom_concat('+', Zero, No_change),
append(Init, [No_change], Temp_list)
)
atom_number(Stoech, Stoech_number)
;
(
findall(N, member([Counter, N], Changes), N_list),
summation(N_list, 1, '', Sum),
append(Init, [Sum], Temp_list)
)
atom_number(Stoech_wp, Stoech_number),
atom_concat('+',Stoech_wp,Stoech)
),
Next_counter is (Counter + 1),
build_change_list(Changes, Next_counter, Length, Temp_list, Result)
)
;
(
Result = Init
).
summation(Input_list, Count, Init, Final) :-
length(Input_list, Length),
nth1(Count, Input_list, To_add),
atom_concat(Init, To_add, Temp),
(
Count < Length
->
(
Next_count is Count + 1,
summation(Input_list, Next_count, Temp, Final)
)
;
(
Final = Temp
)
).
build_stoechiometry(Reactant_list, Product_list, TailSp, TailStoech).
% build_reaction_data(+List_of_reactions, +Species_list, -New_reaction_list)
%
% format the reaction data to the format: [Stoechiometry, Rate]
build_reaction_data([], _Species_list, []) :- !.
build_reaction_data([[RList, PList, Coef]|ReacTail], Species_list, [[Stoech, Coef]|DataTail]) :-
build_stoechiometry(RList, PList, Species_list, Stoech),
build_reaction_data(ReacTail, Species_list, DataTail).
% print_event()
%
% write the events in the output file
print_event(Reactions_data, Stoch_species_order, _Stoch_reactants_order, Hybrid_reactants_list, Constraints_list, Stream) :-
get_option(stochastic_conversion,Rate),
......@@ -1232,10 +1187,10 @@ hybrid_dynamic_simulation(InputFile,OutFilename,Volume,Time,PropTresh,PartTresh)
build_constraints_by_rules(Reactants_by_rules, SortedReactants, 1, [], Hybrid_reactants_list),
build_constraints_list(Hybrid_reactants_list, Constraints_list),
length(SortedSpecies,Species_number),
build_reactants(SortedSpecies, SortedSpecies, Species_number, 1, [], StochReactants),
build_reactants(SortedSpecies, SortedSpecies, StochReactants),
prepare_change_list(Reactions,SortedSpecies,List),
prepare_change_list(Reactions, SortedSpecies, Reaction_list),
build_reaction_data(Reaction_list, SortedSpecies, Reactions_data),
% Because inf/0 is not defined in old versions of SWI-Prolog, we use
% something else
current_prolog_flag(max_tagged_integer, Inf),
......@@ -1243,36 +1198,10 @@ hybrid_dynamic_simulation(InputFile,OutFilename,Volume,Time,PropTresh,PartTresh)
assertz(maximum_reactant_change(0)),
assertz(minimum_product_change(Inf)),
assertz(maximum_product_change(0)),
findall([Species_change, Alpha],
(
member(Element,List),
% To get [RSpecies, PSpecies] from [Reactants, Products]
Element = [Reactants,Products,Alpha],
findall([Rname, Minus_number],
(
member(R, Reactants),
R = Minus_number * Rname
),
RSpecies),
findall([Pname, Add_number],
(
member(P, Products),
P = Add_number * Pname
),
PSpecies),
%== Build list (Species_change) to record species change ==%
particle_sum_list(RSpecies,RSum),
particle_sum_list(PSpecies,PSum),
compute_partial_maximum_change(RSum,PSum),
append_reactants(RSpecies, 1, [], Temp_list, SortedSpecies),
append_products(PSpecies, 1, Temp_list, Changes, SortedSpecies),
build_change_list(Changes, 1, Species_number, [], Species_change)
),
Reactions_data
),
findall([RSpecies,Alpha],
(
member(Element,List),
member(Element,Reaction_list),
Element = [Reactants,_,Alpha],
findall(Rname1,
(
......
:- use_module(hybrid).
:- begin_tests(hybrid).
test(build_constraints_list) :-
hybrid:build_constraints_list([[a],[0],[b,c]], [lower_level1, '0', lower_level3]).
test(build_reaction_data) :-
Sorted_species = [gen,struc,tem,virus],
List = [[[1*gen],[1*tem],c1*'gen_stoch'],[[1*tem],[],c2*('tem_total'/(100*1))],[[1*tem],[1*tem,1*gen],c3*('tem_total'/(100*1))],[[1*gen,1*struc],[1*virus],c4*('gen_stoch'*('struc_total'/(100*1)))]],
hybrid:build_reaction_data(List, Sorted_species, Data),
Data = [[['-1','+0','+1','+0'],c1*'gen_stoch'],[['+0','+0','-1','+0'],c2*('tem_total'/(100*1))],[['+1','+0','+0','+0'],c3*('tem_total'/(100*1))],[['-1','-1','+0','+1'],c4*('gen_stoch'*('struc_total'/(100*1)))]].
:- end_tests(hybrid).
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