1:- module(matrix, [matrix_calc/3, matrix/2]). 2:- use_module(library(lists)). 3:- use_module(util(html)). 4:- use_module(util(misc)). 5:- use_module(util('meta2')). 6:- use_module(util('gb-lex-old')). 7:- use_module(pac(basic)). 8:- use_module(pac(reduce)). 9:- use_module(pac('expand-pac')). 10term_expansion --> pac:expand_pac.
11
12
13:- use_module(pac(op)). 14:- use_module(pac(basic)). 15
16equation_separator(;).
17
24
26
27algebra(X, Y, Z) --> matrix(X, Y, Z).
29is_matrix([_|_]).
30
32subst_rdiv(X, Y) :- subst_functor([((rdiv/2), (/))], X, Y).
33
36
37subst_functor(_, X, X):- (number(X); var(X)), !.
38subst_functor(S, X, Y):- X =..[F|A],
39 maplist(subst_functor(S), A, B),
40 functor(X, F, N),
41 ( memberchk((F/N, G), S)
42 -> H = G
43 ; H = F
44 ),
45 Y =..[H|B].
46
74
75:- bekind(matrix, []). 76square(L) = :vector_to_square@L.
77diag(M) = :diagonal@M.
78unit(N) = :unit_matrix@N0 :- N0 is N.
79inverse(M) = xargs([M0]>>M0^(-1))@M.
80vector(N, A)= V :- length(V, N), maplist(=(A), V).
81t(A) = :transpose_matrix@A.
82tr(A) = :matrix_trace@A.
83-A = scalar(-1)@A.
84A - B = A + scalar(-1, B).
85A + B = :add_m_m@ A@B.
86A * B = :mul_m_m@ A@B.
87s_add(A,B) = :add_s_m(A0)@B :- A0 is A.
88s_mull(A,B) = :mul_s_m(A0)@B :- A0 is A.
89A^N = :exp_m(N0)@ A :- N0 is N.
90max(A, B) = `C :- ( A< B -> C=B; C=A).
91min(A, B) = `C :- ( A< B -> C=A; C=B).
92scalar(A, X)= :mul_s_m(A0)@ X :- A0 is A.
93subst(X, P) = :poly_subst@X@ (math:poly(P)).
94X = `X.
95:- ekind.
96
97% ?- matrix_calc([[[1,2,1],[1,3,1],[-2,-1,0]],[[1,2,1],[1,3,1],[-2,-1,0]],[[1,0,0],[0,1,0],[0,0,1]]], x+y, R).
98% ?- matrix_calc([[[1,2,1],[1,3,1],[-2,-1,0]],[[1,2,1],[1,3,1],[-2,-1,0]],[[1,0,0],[0,1,0],[0,0,1]]], x(1)+x(2), R).
99
100% ?- listing(matrix_calc).
101
102:- bekind(matrix_calc//1,[]).
103x = with([A|_], `A).
104y = with([_,A|_], `A).
105z = with([_,_,A|_], `A).
106x(K) = with(L, :nth1(K, L)).
107square(L) = :vector_to_square@L.
108diag(M) = :diagonal@M.
109unit(N) = :unit_matrix@N0 :- N0 is N.
110inverse(M) = xargs([M0] >> M0^(-1))@M.
111vector(N, A)= V :- length(V, N), maplist(=(A), V).
112t(A) = :transpose_matrix@A.
113tr(A) = :matrix_trace@A.
114-A = scalar(-1)@A.
115A - B = A + scalar(-1, B).
116A + B = :add_m_m@ A@B.
117A * B = :mul_m_m@ A@B.
118s_add(A,B) = :add_s_m(A0)@B :- A0 is A.
119s_mull(A,B) = :mul_s_m(A0)@B :- A0 is A.
120A^N = :exp_m(N0)@ A :- N0 is N.
121max(A, B) = `C :- ( A< B -> C=B; C=A).
122min(A, B) = `C :- ( A< B -> C=A; C=B).
123scalar(A, X)= :mul_s_m(A0)@ X :- A0 is A.
124subst(X, P) = :poly_subst@ X @ (math:poly(P)).
125X