Commit 034cbcf0 authored by Thierry Martinez's avatar Thierry Martinez

Options for numerical simulation

parent ed455f44
......@@ -559,6 +559,19 @@ write_argument(Argument, Doc) :-
<a href="#~a">~a</a><sub><var>n</var></sub>',
[NameId, Name, ValueId, Value, NameId, Name, ValueId, Value]
)
;
Type = '*'(Name: Value)
->
make_id(Name, NameId),
make_id(Value, ValueId),
format(
Doc,
'<a href="#~a">~a</a><sub>1</sub>: <a href="#~a">~a</a><sub>1</sub>,
...,
<a href="#~a">~a</a><sub><var>n</var></sub>:
<a href="#~a">~a</a><sub><var>n</var></sub>',
[NameId, Name, ValueId, Value, NameId, Name, ValueId, Value]
)
;
Type = '*'([Name])
->
......
......@@ -54,7 +54,8 @@ write_headers(Options) :-
write_header('VARIABLE_COUNT', '~d', VariableCount),
write_header('PARAMETER_COUNT', '~d', ParameterCount),
write_option('METHOD', '~a', method, Options),
write_option('ERROR_EPSILON', '~f', error_epsilon, Options),
write_option('ERROR_EPSILON_ABSOLUTE', '~f', error_epsilon_absolute, Options),
write_option('ERROR_EPSILON_RELATIVE', '~f', error_epsilon_relative, Options),
write_option('INITIAL_STEP_SIZE', '~f', initial_step_size, Options),
write_option('TIME_INITIAL', '~f', time_initial, Options),
write_option('TIME_FINAL', '~f', time_final, Options).
......
......@@ -2,26 +2,30 @@
:- begin_tests(gsl).
test('van_der_pol', [
true((
Table = [FirstRow | OtherRows],
append(_, [LastRow], OtherRows),
FirstRow = row(FirstTimeStamp, _, _),
LastRow = row(LastTimeStamp, _, _),
FirstTimeStamp < 1e-5,
LastTimeStamp == 100.0
))]) :-
test('van_der_pol') :-
Options = [
equations: [[1], -[0] + p(0) * [1] * (1 - [0] ^ 2)],
initial_values: [1.0, 0.0],
initial_parameter_values: [10],
method: gsl_odeiv2_step_rk8pd,
error_epsilon: 1e-6,
error_epsilon_absolute: 1e-6,
error_epsilon_relative: 1e-6,
initial_step_size: 1e-6,
precision: 5,
time_initial: 0,
time_final: 100
],
solve(Options, Table).
solve(Options, Table),
(
append(_, [LastRow], Table),
LastRow = row(LastTimeStamp, _, _),
LastTimeStamp == 100.0
->
true
;
print(Table),
nl,
fail
).
:- end_tests(gsl).
......@@ -14,26 +14,26 @@ main(void)
double parameters[PARAMETER_COUNT];
const gsl_odeiv2_step_type *T = METHOD;
gsl_odeiv2_step *s = gsl_odeiv2_step_alloc(T, VARIABLE_COUNT);
gsl_odeiv2_control *c = gsl_odeiv2_control_y_new(ERROR_EPSILON, 0.0);
gsl_odeiv2_evolve * e = gsl_odeiv2_evolve_alloc(VARIABLE_COUNT);
gsl_odeiv2_system sys = {functions, jacobian, VARIABLE_COUNT, parameters};
gsl_odeiv2_driver *d = gsl_odeiv2_driver_alloc_y_new(
&sys, METHOD, INITIAL_STEP_SIZE, ERROR_EPSILON_ABSOLUTE,
ERROR_EPSILON_RELATIVE);
double t = TIME_INITIAL;
double h = INITIAL_STEP_SIZE;
double x[VARIABLE_COUNT];
gsl_odeiv2_step_set_driver(s, d);
gsl_odeiv2_driver_set_nmax(d, 1);
initial_values(x);
initial_parameter_values(parameters);
print_headers(csv);
while (t < TIME_FINAL) {
int status = gsl_odeiv2_evolve_apply(
e, c, s, &sys, &t, TIME_FINAL, &h, x);
if (status != GSL_SUCCESS) {
int status = gsl_odeiv2_driver_apply(d, &t, TIME_FINAL, x);
if (status != GSL_SUCCESS && status != GSL_EMAXITER) {
fprintf(stderr, "error, return value=%d\n", status);
return EXIT_FAILURE;
}
print(csv, t, x);
}
gsl_odeiv2_evolve_free(e);
gsl_odeiv2_control_free(c);
gsl_odeiv2_driver_free(d);
gsl_odeiv2_step_free(s);
return EXIT_SUCCESS;
}
......@@ -50,9 +50,10 @@ solve(Time) :-
initial_values: InitialValues,
initial_parameter_values: InitialParameterValues,
method: gsl_odeiv2_step_rk8pd,
error_epsilon: 1e-6,
error_epsilon_absolute: 1e-6,
error_epsilon_relative: 1e-6,
initial_step_size: 1e-6,
precision: 5,
precision: 6,
time_initial: 0,
time_final: Time
],
......
......@@ -24,6 +24,7 @@
** Numerical temporal properties
* Commands at Top-level
- toplevel.pl
- options.pl
** Loading, listing, importing and exporting models
*** Biocham files
- models.pl
......
......@@ -122,6 +122,12 @@ check_type('='(LeftType, RightType), Equal, NewLeft = NewRight) :-
check_type(LeftType, Left, NewLeft),
check_type(RightType, Right, NewRight).
check_type(':'(LeftType, RightType), Equal, (NewLeft: NewRight)) :-
!,
Equal = (Left: Right),
check_type(LeftType, Left, NewLeft),
check_type(RightType, Right, NewRight).
check_type({SubType}, Set, NewList) :-
!,
set_to_list(Set, 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