Commit d0ca9607 authored by Thierry Martinez's avatar Thierry Martinez

FOLTL

parent 630235db
......@@ -22,10 +22,8 @@ add_conservation(Conservation) :-
simulation, one of those variables will be
eliminated thanks to this conservation law. Be careful if you declare
conservation laws and then plot the result of a previous simulation, the
caption might not be correct.'
),
doc(
'When added, the conservation law will be checked against the reactions
caption might not be correct.
When added, the conservation law will be checked against the reactions
(i.e. purely from stoichiometry), if that fails against the kinetics.
Since these checks are not complete, even a failure will be accepted with
a warning.'
......
:- module(
foltl,
[
validity_domain/1,
validity_domain/2,
op(700, xfy, '<>'),
op(800, xfy, '/\\'),
op(900, xfy, '\\/')
]
).
:- devdoc('\\section{Commands}').
validity_domain(Formula) :-
biocham_command,
type(Formula, foltl),
doc('
solves a FOLTL(R) constraint on the current trace, i.e. computes
the validity domain for the free variables
that make the formula true on the numerical trace.
'),
validity_domain(Formula, Domain),
format('~w\n', [Domain]).
:- devdoc('\\section{Public API}').
validity_domain(Formula, Domain) :-
ExecutableFilename = 'check',
generate_cpp_file(Formula, 'check.inc'),
compile_cpp_program(ExecutableFilename),
export_table('check.csv'),
call_subprocess(
ExecutableFilename, ['check.csv'], [stdout(pipe(ResultStream))]
),
read_term(
ResultStream,
DomainRaw,
[variables(Variables), variable_names(VariableNames)]),
name_variables_and_anonymous(Variables, VariableNames),
!,
reformat_domain(DomainRaw, Domain).
:- devdoc('\\section{Private predicates}').
generate_cpp_file(Formula, Filename) :-
setup_call_cleanup(
open(Filename, write, Stream),
with_output_to(
Stream,
generate_cpp(Formula)
),
close(Stream)
).
:- dynamic(free_variable/2).
:- dynamic(sub_formula/2).
generate_cpp(Formula) :-
set_counter(free_variables, 0),
retractall(free_variable(_, _)),
set_counter(sub_formulae, 0),
retractall(sub_formula(_, _)),
write(
'\c
Domain
compute_domain(Table table) {
'),
declare_free_variables(Formula),
declare_formula(Formula),
write(
' \c
for (Table::reverse_iterator i = table.rbegin(); i != table.rend(); ++i) {
'),
\+ (
sub_formula(SubFormula, Index),
\+ (
generate_sub_formula(SubFormula, Index)
)
),
\+ (
sub_formula(_SubFormula, Index),
\+ (
format(' domain~d = next_domain~d;\n', [Index, Index])
)
),
write(
' \c
}
return domain0;
}
').
generate_sub_formula(Predicate, Index) :-
predicate_binary(Predicate, A, Op, B),
!,
with_output_to(
atom(Left),
generate_expression(A)
),
with_output_to(
atom(Right),
generate_expression(B)
),
format(atom(Constraint), '~a ~a ~a', [Left, Op, Right]),
peek_count(free_variables, FreeVariableCount),
(
(
has_variable(A)
;
has_variable(B)
)
->
format(
' \c
Domain next_domain~d(~d, PPL::EMPTY);
PPL::Constraint_System cs;
cs.insert(~a);
next_domain~d.add_disjunct(PPL::C_Polyhedron(cs));
',
[Index, FreeVariableCount, Constraint, Index]
)
;
format(
' \c
Domain next_domain~d(~d, ~a ? PPL::UNIVERSE : PPL::EMPTY);
',
[Index, FreeVariableCount, Constraint]
)
).
generate_sub_formula('X'(A), Index) :-
sub_formula(A, AIndex),
format(
' \c
Domain next_domain~d(domain~d);
',
[Index, AIndex]
).
generate_sub_formula('F'(A), Index) :-
sub_formula(A, AIndex),
format(
' \c
Domain next_domain~d(domain~d);
for (
Domain::iterator j = next_domain~d.begin();
j != next_domain~d.end(); ++j
) {
next_domain~d.add_disjunct(j->pointset());
}
',
[Index, Index, AIndex, AIndex, Index]
).
generate_sub_formula('G'(A), Index) :-
sub_formula(A, AIndex),
format(
' \c
Domain next_domain~d(domain~d);
next_domain~d.intersection_assign(next_domain~d);
',
[Index, Index, AIndex, AIndex, Index]
).
generate_sub_formula('U'(A, B), Index) :-
generate_sub_formula_u_or_w(A, B, Index).
generate_sub_formula('W'(A, B), Index) :-
generate_sub_formula_u_or_w(A, B, Index).
generate_sub_formula(A /\ B, Index) :-
sub_formula(A, AIndex),
sub_formula(B, BIndex),
format(
' \c
Domain next_domain~d(next_domain~d);
next_domain~d.intersection_assign(next_domain~d);
',
[Index, AIndex, Index, BIndex]
).
generate_sub_formula(A \/ B, Index) :-
sub_formula(A, AIndex),
sub_formula(B, BIndex),
format(
' \c
Domain next_domain~d(next_domain~d);
for (
Domain::iterator j = next_domain~d.begin();
j != next_domain~d.end(); ++j
) {
next_domain~d.add_disjunct(j->pointset());
}
',
[Index, AIndex, BIndex, BIndex, Index]
).
generate_sub_formula_u_or_w(A, B, Index) :-
sub_formula(A, AIndex),
sub_formula(B, BIndex),
format(
' \c
Domain next_domain~d(domain~d);
next_domain~d.intersect_assign(next_domain~d);
for (
Domain::iterator j = next_domain~d.begin();
j != next_domain~d.end(); ++j
) {
next_domain~d.add_disjunct(j->pointset());
}
',
[Index, Index, Index, AIndex, BIndex, BIndex, Index]
).
generate_expression(Number) :-
number(Number),
!,
format('~w', [Number]).
generate_expression(A) :-
free_variable(A, Index),
!,
format('x~d', [Index]).
generate_expression(A) :-
get_current_table(Table),
columns(Table, Index, A),
!,
format('(*i)[~d]', [Index]).
has_variable(A) :-
free_variable(A, _Index).
declare_free_variables(Predicate) :-
predicate_binary(Predicate, A, _Op, B),
!,
declare_free_variables_expression(A),
declare_free_variables_expression(B).
declare_free_variables(Formula) :-
foltl_unary(Formula, SubFormula),
!,
declare_free_variables(SubFormula).
declare_free_variables(Formula) :-
foltl_binary(Formula, Formula0, Formula1),
!,
declare_free_variables(Formula0),
declare_free_variables(Formula1).
declare_free_variables_expression(X) :-
number(X),
!.
declare_free_variables_expression(X) :-
atom(X),
!,
(
get_current_table(Table),
columns(Table, _Index, X)
->
true
;
free_variable(X, _Index)
->
true
;
count(free_variables, Index),
asserta(free_variable(X, Index)),
format(' PPL::Variable x~d(~d);\n', [Index, Index])
).
declare_formula(Formula) :-
(
(
Formula = 'G'(_)
;
Formula = 'W'(_, _)
)
->
Initial = 'PPL::UNIVERSE'
;
Initial = 'PPL::EMPTY'
),
index_formula(Formula, Initial),
declare_sub_formulae(Formula).
declare_sub_formulae(Formula) :-
foltl_atomic(Formula),
!.
declare_sub_formulae(Formula) :-
foltl_unary(Formula, SubFormula),
!,
declare_formula(SubFormula).
declare_sub_formulae(Formula) :-
foltl_binary(Formula, Formula0, Formula1),
!,
declare_formula(Formula0),
declare_formula(Formula1).
predicate_binary(A = B, A, '==', B).
predicate_binary(A <> B, A, '!=', B).
predicate_binary(A < B, A, '<', B).
predicate_binary(A <= B, A, '<=', B).
predicate_binary(A > B, A, '>', B).
predicate_binary(A >= B, A, '>=', B).
foltl_atomic(A) :-
predicate_binary(A, _, _, _).
foltl_unary('X'(Formula), Formula).
foltl_unary('F'(Formula), Formula).
foltl_unary('G'(Formula), Formula).
foltl_binary('U'(Formula0, Formula1), Formula0, Formula1).
foltl_binary('W'(Formula0, Formula1), Formula0, Formula1).
foltl_binary(Formula0 /\ Formula1, Formula0, Formula1).
foltl_binary(Formula0 \/ Formula1, Formula0, Formula1).
index_formula(Formula, Initial) :-
(
sub_formula(Formula, _Index)
->
true
;
peek_count(free_variables, FreeVariableCount),
count(sub_formulae, Index),
asserta(sub_formula(Formula, Index)),
format(
' Domain domain~d(~d, ~a);\n',
[Index, FreeVariableCount, Initial]
)
).
compile_cpp_program(ExecutableFilename) :-
library_path(LibraryPath),
absolute_file_name(
'ppl_validity_domain.cc', ValidityDomainCc, [relative_to(LibraryPath)]
),
call_subprocess(
path('g++'),
['-o', ExecutableFilename, '-I.', ValidityDomainCc,
'-lppl_c', '-lppl', '-lgmpxx', '-lgmp'],
[]).
reformat_domain({ Conjunction0 }, Conjunction1) :-
reformat_domain_conjunction(Conjunction0, Conjunction1).
reformat_domain((A0, B0), A1 \/ B1) :-
reformat_domain(A0, A1),
reformat_domain(B0, B1).
reformat_domain(false, false).
reformat_domain(true, true).
reformat_domain_conjunction(A0 = B0, A1 = B1) :-
reformat_expression(A0, A1),
reformat_expression(B0, B1).
reformat_domain_conjunction(A0 < B0, A1 < B1) :-
reformat_expression(A0, A1),
reformat_expression(B0, B1).
reformat_domain_conjunction(A0 <= B0, A1 <= B1) :-
reformat_expression(A0, A1),
reformat_expression(B0, B1).
reformat_domain_conjunction(A0 > B0, A1 > B1) :-
reformat_expression(A0, A1),
reformat_expression(B0, B1).
reformat_domain_conjunction(A0 >= B0, A1 >= B1) :-
reformat_expression(A0, A1),
reformat_expression(B0, B1).
reformat_domain_conjunction((A0, B0), A1 /\ B1) :-
reformat_domain_conjunction(A0, A1),
reformat_domain_conjunction(B0, B1).
reformat_expression(Number, Number) :-
number(Number),
!.
reformat_expression(DimensionName, ColumnName) :-
atom(DimensionName),
!,
atom_length(DimensionName, NameLength),
(
NameLength > 1
->
sub_atom(DimensionName, 0, 1, _, Letter),
sub_atom(DimensionName, 1, _, 0, Counter),
letter_index(Letter, LetterIndex),
Index is Counter * 26 + LetterIndex
;
letter_index(DimensionName, Index)
),
get_current_table(Table),
once(columns(Table, Index, ColumnName)).
letter_index(Letter, Index) :-
once(sub_atom('ABCDEFGHIJKLMNOPQRSTUVWXYZ', Index, 1, _, Letter)).
:- use_module(library(plunit)).
:- begin_tests(foltl).
test('validity_domain F', [true(Domain == ('x' >= 10))]) :-
clear_model,
add_table(table, [row('#x'), row(20), row(10), row(30), row(15)]),
validity_domain('F'('x' <= 'v'), Domain).
test(
'validity_domain =',
[true(Domain == (x = 15 \/ x = 30 \/ x = 10 \/ x = 20))]
) :-
clear_model,
add_table(table, [row('#x'), row(20), row(10), row(30), row(15)]),
validity_domain('F'(x = v), Domain).
test('validity_domain X', [true(Domain == (x = 30 \/ x = 10))]) :-
clear_model,
add_table(table, [row('#x'), row(20), row(10), row(30), row(15)]),
validity_domain('F'(x = v /\ 'X'(x > 10)), Domain).
:- end_tests(foltl).
......@@ -258,13 +258,11 @@ list_attributes(GraphObject) :-
get_current_graph(Id) :-
find_item([kind: current_graph, item: Id]).
get_selection(current_graph, [Id]).
set_current_graph(Id) :-
deselect_graph,
add_item([kind: current_graph, item: Id, id: CurrentGraphId]),
add_dependency(CurrentGraphId, Id).
set_selection(current_graph, [Id]).
get_graph_name(Id, Name) :-
......@@ -322,10 +320,6 @@ add_graph_object(GraphObject, Id) :-
).
deselect_graph :-
delete_items([kind: current_graph]).
isolated_vertices(GraphId, VerticesId) :-
findall(
VertexId,
......
......@@ -21,7 +21,7 @@ test('set_graph_name', [true(Graphs == [my_graph])]) :-
set_graph_name(my_graph),
all_items([kind: graph], Graphs).
test('list_graphs', [true(Atom == '[0] new_graph\n')]) :-
test('list_graphs', [true(Atom == '[0]*new_graph\n')]) :-
clear_model,
new_graph,
with_output_to(atom(Atom), list_graphs).
......
......@@ -20,7 +20,7 @@ solve(Options, Table) :-
ExecutableFilename = 'ode',
write_gsl_file(CFilename, Options),
compile_c_program(ExecutableFilename),
call_subprocess(ExecutableFilename, []),
call_subprocess(ExecutableFilename, [], []),
csv_read_file('ode.csv', Table).
......@@ -293,11 +293,6 @@ compile_c_program(ExecutableFilename) :-
absolute_file_name('gsl_solver.c', GslSolverC, [relative_to(LibraryPath)]),
call_subprocess(
path(gcc),
['-o', ExecutableFilename, '-I', '.', GslSolverC,
'-lgsl', '-lgslcblas', '-lm']).
call_subprocess(ExecutableFilename, Arguments) :-
absolute_file_name(ExecutableFilename, Absolute),
process_create(Absolute, Arguments, [process(Pid)]),
process_wait(Pid, exit(0)).
['-o', ExecutableFilename, '-I.', GslSolverC,
'-lgsl', '-lgslcblas', '-lm'],
[]).
#include <iostream>
#include <sstream>
#include <fstream>
#include <string>
#include <stdlib.h>
#include <vector>
#include <ppl.hh>
namespace PPL = Parma_Polyhedra_Library;
typedef std::vector<double> Row;
typedef std::vector<Row> Table;
Table
read_CSV_stream(std::istream &stream)
{
Table result;
std::string line;
while (std::getline(stream, line)) {
if (line[0] == '#') {
continue;
}
std::stringstream lineStream(line);
Row row;
std::string cell;
while (std::getline(lineStream, cell, ',')) {
row.push_back(atof(cell.c_str()));
}
result.push_back(row);
}
return result;
}
typedef PPL::Pointset_Powerset<PPL::C_Polyhedron> Domain;
#include "check.inc"
int
main(int argc, char *argv[])
{
std::ifstream file(argv[1]);
Table table = read_CSV_stream(file);
Domain result = compute_domain(table);
using PPL::IO_Operators::operator<<;
std::cout << result << "." << std::endl;
return EXIT_SUCCESS;
}
......@@ -36,7 +36,9 @@
add_file_suffix/2,
inherits/1,
inherits/2,
begin_command/0
begin_command/0,
get_selection/2,
set_selection/2
]
).
......@@ -334,7 +336,7 @@ load_biocham_stream(Stream) :-
current_models(Models) :-
nb_getval(current_models, Models).
get_selection(current_models, Models).
set_current_models(Models) :-
......@@ -343,7 +345,7 @@ set_current_models(Models) :-
->
throw(error(cannot_select_no_models))
;
nb_setval(current_models, Models)
set_selection(current_models, Models)
).
......@@ -641,7 +643,14 @@ list_ids_aux(Options, Ids) :-
assertz(listed_item(Counter, Id)),
item(Id, _, _, Item),
indent(Indent),
format('[~d] ~w\n', [Counter, Item]),
(
selection(_, Id)
->
Selected = '*'
;
Selected = ' '
),
format('[~d]~a~w\n', [Counter, Selected, Item]),
(
Recursive = yes
->
......@@ -692,23 +701,16 @@ delete_item(Id) :-
)
)
),
retractall(selection(_, Id)),
(
KindItem = model
->
current_models(CurrentModels),
(
select(Id, CurrentModels, OtherModels)
selection(current_models, _)
->
nb_setval(current_models, OtherModels),
(
OtherModels = []
->
new_model
;
true
)
;
true
;
new_model
)
;
true
......@@ -816,3 +818,29 @@ inherits_from(Id, AncestorId) :-
inherits_from(Id, AncestorId) :-
directly_inherits_from(Id, IntermediateAncestorId),
inherits_from(IntermediateAncestorId, AncestorId).
:- dynamic(selection/2).
get_selection(SelectionName, Ids) :-
findall(
Id,
selection(SelectionName, Id),
Ids).
set_selection(SelectionName, IdsOrOptions) :-
retractall(selection(SelectionName, _)),
\+ (
(
IdsOrOptions = [_: _| _]
->
item([id: Id | IdsOrOptions])
;
member(Id, IdsOrOptions)
),
\+ (
assertz(selection(SelectionName, Id))
)
).
......@@ -56,7 +56,7 @@ solve(Time) :-
time_final: Time
],
solve(Options, Table),
add_trace('numerical_simulation', Table).
add_table('numerical_simulation', Table).
gather_headers(Headers) :-
......
......@@ -61,8 +61,8 @@ export_plot_to_png(OutputFile) :-
export_plot(FileTemplate, Options) :-
format(atom(PlotFile), '~a.plot', [FileTemplate]),
format(atom(CsvFile), '~a.csv', [FileTemplate]),
export_trace(CsvFile),
get_trace_data([HeaderRow | _ ]),
export_table(CsvFile),
get_table_data([HeaderRow | _ ]),
HeaderRow =.. [row, _TimeHeader | Headers],
setup_call_cleanup(
open(PlotFile, write, Stream),
......
:- module(
tables,
[
load_table/1,
export_table/1,
select_table/1,
rename_table/1,
list_tables/0,
delete_table/1,
list_columns/0,
column/1,
delete_column/1,
rename_column/2,
add_table/2,
get_table_data/1,
get_current_table/1,
columns/3
]
).
:- devdoc('\\section{Commands}').
load_table(InputFile) :-