Commit 34d482bc authored by HEMERY Mathieu's avatar HEMERY Mathieu
Browse files

Rewrite various print_dynamic_event predicates

parent eb012e7a
......@@ -321,23 +321,6 @@ build_constraints_list([List| InTail], Counter, [Constraint| OutTail]) :-
),
Counterp is Counter + 1,
build_constraints_list(InTail, Counterp, OutTail).
% build_reactants(Species, ODESpecies, Reactants)
%
% Construct the list of reactant species in an adequate format
build_reactants([], _ODESpecies, []).
build_reactants([Species|SpTail], ODESpecies, [Reactant|ReacTail]) :-
(
member(Species, ODESpecies)
->
species_to_total(Species, Reactant)
;
species_to_stoch(Species, Reactant)
),
build_reactants(SpTail, ODESpecies, ReacTail).
%! write_hybrid_ode(+Stream, +Hybrid_list)
......@@ -893,7 +876,6 @@ 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),
build_reactants(SortedSpecies, SortedSpecies, StochReactants),
prepare_change_list(Reactions, SortedSpecies, Reaction_list),
build_reaction_data(Reaction_list, SortedSpecies, Reactions_data),
......@@ -901,7 +883,7 @@ hybrid_dynamic_simulation(InputFile,OutFilename,Volume,Time,PropTresh,PartTresh)
compute_maximum_change(Reaction_list, MaxParticleChange),
ParticuleTreshold is MaxParticleChange*PartTresh,
print_dynamic_event(SortedSpecies,Reactions_data, HybridDynamicReactants,ParticuleTreshold,PropTresh, SortedStochSpecies, StochReactants,Hybrid_reactants_list, Constraints_list, Stream),
print_dynamic_event(SortedSpecies,Reactions_data, HybridDynamicReactants,ParticuleTreshold,PropTresh, SortedStochSpecies, Hybrid_reactants_list, Constraints_list, Stream),
retractall(reaction_kinetics(_,_)),
format(Stream, "~n% Plot part~n~n", []),
......@@ -966,59 +948,44 @@ build_hybrid_subroutine([_N*X|InputTail], [Name|OutputTail]) :-
build_hybrid_subroutine(InputTail, OutputTail).
print_dynamic_event_condition(SortedSpecies,HybridDynamicReactants,ReactionNumber,PropTresh,MaxParticleChange,Stream) :-
nth1(ReactionNumber,HybridDynamicReactants,ReactionInfo),
nth1(1,ReactionInfo,Reactants),
write(Stream,' k'),
write(Stream,ReactionNumber),
write(Stream,'_diff = (if '),
nb_getval(reactions,Reactions),
nth1(ReactionNumber,Reactions,Reaction),
convert_MA2(Reaction,SortedSpecies,NewKinetics,_,_),
write(Stream,NewKinetics),
write(Stream, ' >= 1/(step_size * '),
write(Stream,PropTresh),
write(Stream,')'),
print_dynamic_event_condition_ifThenElse(Reactants,MaxParticleChange,Stream),
write(Stream,' then 1 else 0),'),
nl(Stream),
(
length(HybridDynamicReactants,NumberOfReactions),
ReactionNumber < NumberOfReactions
->
NextReactionNumber is ReactionNumber + 1,
print_dynamic_event_condition(SortedSpecies,HybridDynamicReactants,NextReactionNumber,PropTresh,MaxParticleChange,Stream)
;
true
).
print_dynamic_event_condition_ifThenElse([Reactant|Reactants],MaxParticleChange,Stream) :-
write(Stream,' and '),
write(Stream,Reactant),
write(Stream,' >= '),
write(Stream,MaxParticleChange),
print_dynamic_event_condition_ifThenElse(Reactants,MaxParticleChange,Stream).
print_dynamic_event_condition_ifThenElse([],_,_).
%! print_dynamic_event()
%
%
print_dynamic_event(SortedSpecies,Reactions_data, HybridDynamicReactants,MaxParticleChange,PropTresh,Stoch_species_order, _Stoch_reactants_order, Hybrid_reactants_list, Constraints_list, Stream) :-
print_dynamic_event(SortedSpecies, Reactions_data, HybridDynamicReactants,MaxParticleChange,PropTresh,Stoch_species_order, Hybrid_reactants_list, Constraints_list, Stream) :-
get_option(stochastic_conversion,Rate),
length(Reactions_data, Number_of_alpha),
print_alpha_parameters(Number_of_alpha, Stream),
print_constraints_parameters(Constraints_list, Stream),
% 1. Build Starting String for the event.
write(Stream, 'add_event(Time > tau,'),nl(Stream),
print_dynamic_event_condition(SortedSpecies,HybridDynamicReactants,1,PropTresh,MaxParticleChange,Stream),
format(Stream, 'add_event(Time > tau,~n', []),
print_dynamic_event_condition(SortedSpecies,HybridDynamicReactants,PropTresh,MaxParticleChange,Stream),
print_static_alpha(Stream),
write(Stream,' delta_t = (if alpha_sum <= 0 then step_size else (-1 / alpha_sum) * log(random_float)),'),nl(Stream),
write(Stream,' tau = Time+delta_t/'),write(Stream,Rate),write(Stream,','),nl(Stream),
format(Stream, " tau = Time+delta_t/~w,~n", [Rate]),
write(Stream,' ran = random_float,'),nl(Stream),
print_hybrid_species_constraint(Hybrid_reactants_list, 1, Stream),
build_alphas(Number_of_alpha, Alpha_accumulation_list),
% 3. Print reactions in event
% Build the print_helper_list.
build_normal(Reactions_data, Print_helper_list, Alpha_accumulation_list, Constraints_list, Stoch_species_order, 1),
print_change_list(Print_helper_list, Stoch_species_order, Stream).
print_dynamic_event_condition(Sp_list, Reac_info_list, PropTresh, Max_change, Stream) :-
nb_getval(reactions, Reac_list),
print_dynamic_event_condition(Sp_list, Reac_info_list, Reac_list, 1, PropTresh, Max_change, Stream).
print_dynamic_event_condition(_Species, [], [], _Counter, _PropTresh, _Max_change, _Stream).
print_dynamic_event_condition(Species, [[Reactants|_]|Reac_info_tail], [Reac|Reac_tail], Counter,PropTresh,Max_change,Stream) :-
format(Stream, " k~d_diff = (if ", [Counter]),
convert_MA2(Reac, Species, NewKinetics, _, _),
format(Stream, "~w >= 1/(step_size * ~w)", [NewKinetics, PropTresh]),
print_dynamic_event_condition_ifThenElse(Reactants,Max_change,Stream),
format(Stream,' then 1 else 0),~n', []),
NextCounter is Counter + 1,
print_dynamic_event_condition(Species,Reac_info_tail, Reac_tail, NextCounter,PropTresh,Max_change,Stream).
print_dynamic_event_condition_ifThenElse([],_,_).
print_dynamic_event_condition_ifThenElse([Reactant|Tail],Max_change,Stream) :-
format(Stream, ' and ~w >= ~w', [Reactant, Max_change]),
print_dynamic_event_condition_ifThenElse(Tail,Max_change,Stream).
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