Skip to content

GitLab

Commit 9cdf8bdc authored by Mathieu Hemery's avatar Mathieu Hemery
Browse files

Rewrite eval_coeff using substitute tool

parent a2643f19
Hide whitespace changes
Inline Side-by-side
rosenbrock.pl
View file @ 9cdf8bdc
......@@ -4,6 +4,11 @@
]
).
% Only for separate compilation/linting
:- use_module(doc).
:- use_module(arithmetic_rules). % simplify,
:- use_module(util). % substitute,
%! rosenbrock_init
%
% Initialize coefficients for the rosenbrock method, the errors and the step size
......@@ -122,8 +127,8 @@ ode_solver([
% Sublist is the N first element of List
take(Size, List, Sublist) :-
prefix(Sublist, List),
length(Sublist, Size).
length(Sublist, Size),
prefix(Sublist, List).
%! initialize_functions_list
......@@ -195,10 +200,9 @@ initialize_functions :-
append(VariablesNames,FunctionsVariables,NewVariables),
findall(
InitialValue,
(member(function(_,Expression),Functions),
eval_coeff(Expression,InitialState,Parameters,0,InitialValue1),
simplify(InitialValue1,InitialValue2),
InitialValue is InitialValue2
(
member(function(_,Expression),Functions),
eval_coeff(Expression,InitialState,Parameters,0,InitialValue)
),
InitialValues
),
......@@ -222,9 +226,7 @@ update_functions_state(CurrentState,Time) :-
FunctionValue,
(
member(function(_,Expression),Functions),
eval_coeff(Expression,CurrentState,Parameters,Time,FunctionValue1),
simplify(FunctionValue1,FunctionValue2),
FunctionValue is FunctionValue2
eval_coeff(Expression,CurrentState,Parameters,Time,FunctionValue)
),
FunctionsValues
),
......@@ -232,97 +234,52 @@ update_functions_state(CurrentState,Time) :-
NewLastRow =.. [row,Time|NewCurrentState],
append(Table,[NewLastRow],NewTable),
nb_setval(numerical_table,NewTable).
/* Evaluates a matrix/row/coefficient with informations on current state of the simulation and the current parameters */
eval_matrix([],_,_,[]).
eval_matrix([Row|T],CurrentState,Parameters,[RowEvaluated|ET]) :- % CurrentState is a list like [1,5.26,-2,0], the concentrations of the variables [0], [1], [2] etc...
eval_row(Row,CurrentState,Parameters,RowEvaluated),
eval_matrix(T,CurrentState,Parameters,ET).
eval_row([],_,_,[]).
eval_row([Element|T],CurrentState,Parameters,[ElementEvaluated2|ET]) :-
eval_coeff(Element,CurrentState,Parameters,_,ElementEvaluatedUnSimplified),!,
simplify(ElementEvaluatedUnSimplified,ElementEvaluated1),
ElementEvaluated2 is ElementEvaluated1, % There are issues with - (a*b), so it prevents these kinds of expression to appear in the Jacobian.
eval_row(T,CurrentState,Parameters,ET).
eval_coeff(Element1 ^ Element2,CurrentState,Parameters,Time,ElementEvaluated1 ^ ElementEvaluated2) :-
eval_coeff(Element1,CurrentState,Parameters,Time,ElementEvaluated1),
eval_coeff(Element2,CurrentState,Parameters,Time,ElementEvaluated2).
eval_coeff(Element1 + Element2,CurrentState,Parameters,Time,ElementEvaluated1 + ElementEvaluated2) :-
eval_coeff(Element1,CurrentState,Parameters,Time,ElementEvaluated1),
eval_coeff(Element2,CurrentState,Parameters,Time,ElementEvaluated2).
eval_coeff(Element1 - Element2,CurrentState,Parameters,Time,ElementEvaluated1 - ElementEvaluated2) :-
eval_coeff(Element1,CurrentState,Parameters,Time,ElementEvaluated1),
eval_coeff(Element2,CurrentState,Parameters,Time,ElementEvaluated2).
eval_coeff(Element1 * Element2,CurrentState,Parameters,Time,ElementEvaluated1 * ElementEvaluated2) :-
eval_coeff(Element1,CurrentState,Parameters,Time,ElementEvaluated1),
eval_coeff(Element2,CurrentState,Parameters,Time,ElementEvaluated2).
eval_coeff(Element1 / Element2,CurrentState,Parameters,Time,ElementEvaluated1 / ElementEvaluated2) :-
eval_coeff(Element1,CurrentState,Parameters,Time,ElementEvaluated1),
eval_coeff(Element2,CurrentState,Parameters,Time,ElementEvaluated2).
eval_coeff(-Element,CurrentState,Parameters,Time,-ElementEvaluated) :-
eval_coeff(Element,CurrentState,Parameters,Time,ElementEvaluated).
eval_coeff([IndexVariable],CurrentState,_,_,ElementEvaluated) :-
nth0(IndexVariable,CurrentState,ElementEvaluated).
eval_coeff(Element,_,_,_,Element) :-
number(Element).
eval_coeff(p(ParameterIndex),CurrentState,Parameters,Time,ElementEvaluated) :-
nth0(ParameterIndex,Parameters,Element2),
eval_coeff(Element2,CurrentState,Parameters,Time,ElementEvaluated).
eval_coeff(floor(Element),CurrentState,Parameters,Time,Result) :-
!,
eval_coeff(Element,CurrentState,Parameters,Time,ElementEval),
Result is floor(ElementEval).
eval_coeff(random,_,_,_,Result) :-
!,Result is random_float.
eval_coeff(infinity,_,_,_,Result) :-
!,Result is inf.
eval_coeff(t,_,_,Time,Time) :-
!.
eval_coeff(step_size,_,_,_,Result) :-
!,
nb_getval(hdid,Result).
eval_coeff(min(Element1, Element2),CurrentState,Parameters,Time,min(ElementEvaluated1,ElementEvaluated2)) :-
!,
eval_coeff(Element1,CurrentState,Parameters,Time,ElementEvaluated1),
eval_coeff(Element2,CurrentState,Parameters,Time,ElementEvaluated2).
eval_coeff(max(Element1, Element2),CurrentState,Parameters,Time,max(ElementEvaluated1,ElementEvaluated2)) :-
!,
eval_coeff(Element1,CurrentState,Parameters,Time,ElementEvaluated1),
eval_coeff(Element2,CurrentState,Parameters,Time,ElementEvaluated2).
eval_coeff(log(Element),CurrentState,Parameters,Time,log(ElementEvaluated)) :-
!,
eval_coeff(Element,CurrentState,Parameters,Time,ElementEvaluated).
%! eval_matrix(+Symbolic_Matrix, +Current_State, +Parameters, -Evaluated_Matrix)
%
% Evaluates matrix with informations on current state of the simulation and parameters
% Current_State is the list of the species concentration at current time
eval_coeff(exp(Element),CurrentState,Parameters,Time,exp(ElementEvaluated)) :-
!,
eval_coeff(Element,CurrentState,Parameters,Time,ElementEvaluated).
eval_matrix([],_,_,[]).
eval_matrix([Row|T],CurrentState,Parameters,[RowEvaluated|ET]) :-
eval_row(Row,CurrentState,Parameters,RowEvaluated),
eval_matrix(T,CurrentState,Parameters,ET).
eval_coeff(cos(Element),CurrentState,Parameters,Time,cos(ElementEvaluated)) :-
!,
eval_coeff(Element,CurrentState,Parameters,Time,ElementEvaluated).
eval_row([],_,_,[]).
eval_row([Element|T],CurrentState,Parameters,[ElementEvaluated|ET]) :-
eval_coeff(Element,CurrentState,Parameters,_,ElementEvaluated),!,
eval_row(T,CurrentState,Parameters,ET).
eval_coeff(Expr, CurrentState, Parameters, Time, Evaluation) :-
length(CurrentState, NSpecies),
generate_variable_index(NSpecies, IndexSpecies),
length(Parameters, Nparam),
generate_param_index(Nparam, IndexParameters),
nb_getval(hdid, StepSize),
append([IndexSpecies, IndexParameters, [t, random, infinity, step_size]], Index),
append([CurrentState, Parameters, [Time, random_float, inf, StepSize]], Values),
substitute(Index, Values, Expr, Expr_trans),
Evaluation is Expr_trans.
generate_variable_index(N, List) :-
Nm is N-1,
numlist(0,Nm,List_tempo),
put_bracket(List_tempo, List).
put_bracket([], []) :- !.
put_bracket([Head|Tail], [[Head]|BTail]) :-
put_bracket(Tail, BTail).
generate_param_index(N, List) :-
Nm is N-1,
numlist(0,Nm,List_tempo),
put_parenthesis(List_tempo, List).
put_parenthesis([], []) :- !.
put_parenthesis([Head|Tail], [p(Head)|PTail]) :-
put_parenthesis(Tail, PTail).
eval_coeff(sin(Element),CurrentState,Parameters,Time,sin(ElementEvaluated)) :-
!,
eval_coeff(Element,CurrentState,Parameters,Time,ElementEvaluated).
eval_state([],_,[]).
......@@ -918,9 +875,7 @@ apply_changes(State,Time, Parameter = Expression) :-
IndexParameter is Parameter + 1,
nb_getval(parameters_list,Parameters),
eval_coeff(Expression,State,Parameters,Time,Value),
simplify(Value,Value1),
Value2 is Value1,
matrix:replace(Parameters,IndexParameter,Value2,NewParameters),
matrix:replace(Parameters,IndexParameter,Value,NewParameters),
%format("AncientParameters are ~w and NewParameters are ~w ~n",[Parameters,NewParameters]),
nb_setval(parameters_list,NewParameters).
......
rosenbrock.plt
View file @ 9cdf8bdc
......@@ -2,4 +2,19 @@
:- begin_tests(rosenbrock).
test(eva_coeff_1) :-
rosenbrock:eval_coeff(1*[0],[1.5,2.5,3.5],[-1,2], 1, R1), R1 =:= 1.5,
rosenbrock:eval_coeff(2*[0] + p(1),[1.5,2.5,3.5],[-1,2], 1, R2), R2 =:= 5.0,
rosenbrock:eval_coeff(t,[1.5,2.5,3.5],[-1,2], 1.23, R3), R3 =:= 1.23,
rosenbrock:eval_coeff(2*[0] - p(0),[1.5,2.5,3.5],[-1,2], 1, R4), R4 =:=4.0,
rosenbrock:eval_coeff(-p(1),[1.5,2.5,3.5],[-1,2], 1, R5), R5 =:= -2,
rosenbrock:eval_coeff(25/[1],[1.5,2.5,3.5],[-1,2], 1, R6), R6 =:= 10,
rosenbrock:eval_coeff([1]^p(1),[1.5,2.5,3.5],[-1,2], 1, R7), R7 =:= 6.25.
test(eva_coeff_2) :-
rosenbrock:eval_coeff(floor(1*[0]),[1.5,2.5,3.5],[-1,2], 1, R1), R1 =:= 1,
rosenbrock:eval_coeff(min([0], p(1)),[1.5,2.5,3.5],[-1,2], 1, R2), R2 =:= 1.5,
rosenbrock:eval_coeff(infinity,[1.5,2.5,3.5],[-1,2], 1, R3), R3 =:= inf,
rosenbrock:eval_coeff(random,[1.5,2.5,3.5],[-1,2], 1, R4).
:- end_tests(rosenbrock).
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