Calcul numeric / simbolic cu polinoame in Prolog
<![CDATA[
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% Calcul numeric / simbolic cu polinoame
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% Paul Mateescu, grupa 309-Informatica, an III-ZI, 2009.
% Facultatea de Matematica si Informatica, Universitatea Spiru Haret Bucuresti.
% http://fmi.spiruharet.ro/bodorin/aicl.html
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% Prolog, TP 2.0
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
/************************
Calculul sumei, produsului si derivatei polinoamelor.
Solutia propusa implica aducerea unui polinom la o forma canonica de
lista de monoame ordonate strict crescator in ordinea gradului.
Acest lucru se realizeaza prin transformarea unui polinom intr-un arbore binar in
ale carui noduri avem monoamele, astfel incat toate monoamele din sub-arborele stang
sa aiba grad <= monomul din nodul considerat, iar cele din sub-arborele drept sa aiba
gradul strict mai mare.
Arborele binar va fi apoi traversat in inordine, rezultand o sortare a monoamelor in
sens crescator al gradelor. Polinomului rezultat i se aplica reguli de simplificare:
monoamele de acelasi grad se aduna, monoamele de coeficient zero se elimina.
Suma a doua polinoame se face "injectand" monoamele ambelor in acelasi arbore binar si
aducand rezultatul la forma canonica.
Produsul se face uzual, ca suma de produs al monoamelor unui polinom cu celalalt polinom
si aducand rezultatul la forma canonica.
Derivarea se face monom cu monom si adaugand rezultatul la forma canonica.
Am introdus si predicate care afiseaza rezultatul intr-o forma mai lizibila, producand afisari de genul:
P(X) = +1+X
Q(X) = +1-X+X^2
P(X)*Q(X) = +1+X^3
Testarea usoara se face folosind predicatele fara argumente:
testsuma
testprod
testderiv
*************************/
DOMAINS
i=integer
r=real
lr=r*
li=i*
c=real % coefficient
v=string % symbolic variable
p=integer % power
mon=m(real,integer) % monomial
monom=mon*
poly=mon* % polynomial
tree = t(mon, tree, tree)
PREDICATES
/************************
* De aici incepe rezolvarea temei (a se vedea adresa http://fmi.spiruharet.ro/bodorin/labs/ai/091109-13-ai-a3-s1-l005-symbolic-computation-polynomials.html ):
*************************/
epsilon(r) % Precizia sub care se considera ca un coeficient este de fapt zero
zero(r) % Reuseste daca argumentul este zero in acceptiunea de mai sus, esueaza altfel
unu(r) % Reuseste daca argumentul este 1 in acceptiunea de mai sus
sum(poly, poly, poly) % Insumeaza primele doua polinoame, returnand rezultatul in forma canonica
prod(poly, poly, poly) % Produsul a doua polinoame, returnand rezultatul in forma canonica
prodh(poly, poly, poly, poly) % Ajutator pentru prod(poly, poly, poly)
prodmon(mon, poly, poly) % Produsul unui monom cu un polinom, returnand rezultatul in forma canonica
prodmonh(mon, poly, poly) % Ajutator pentru prodmon(mon, poly, poly)
deriv(poly, poly) % Deriveaza un polinom, returnand rezultatul in forma canonica
derivh(poly, poly) % Ajutator pentru deriv(poly, poly)
proper(poly, poly) % Transforma un polinom in forma canonica
lesseq(mon, mon) % Primul argument este <= cu al doilea daca gradele lor sunt in relatia "<="
maketree(poly, tree) % Transforma un polinom intr-un arbore binar al monoamelor componente, sortate in arbore conform relatiei "<="
inject(tree, mon) % Injecteaza un monom in arborele binar, conform relatiei "<="
inorder(tree, poly) % Traverseaza in inordine arborele binar, returnand polinomul in forma canonica
inorderh(tree, poly, poly) % Predicat ajutator pentru inorder(tree, poly)
simpl(poly, poly) % Simplifica un polinom pentru a-l aduce la forma canonica
simplh(poly, poly) % Predicat ajutator pentru simplifica(poly, poly)
writepoly(poly) % Scrie un polinom intr-o forma cat mai lizibila
writemon(mon) % Scrie un monom intr-o forma cat mai lizibila
testsuma % Predicat de test, cere 2 polinoame, le afiseaza in forma canonica si suma lor
testprod % Predicat de test, cere 2 polinoame, le afiseaza in forma canonica si produsul lor
testderiv % Predicat de test, cere un polinom, il afiseaza in forma canonica si derivata sa
CLAUSES
/***************************************************************/
% Precizia sub care se considera ca un numar este egal cu zero
% (o)
epsilon(1E-16).
/***************************************************************/
% Verifica daca numarul este aproximativ egal cu 0, cu toleranta EPS
% (i)
zero(X) :-
epsilon(EPS),
abs(X) <= EPS.
/***************************************************************/
% Verifica daca numarul este aproximativ egal cu 1, cu toleranta EPS
% (i)
unu(X) :-
Z = X - 1.0,
zero(Z).
/***************************************************************/
% Insumeaza doua polinoame:
% (i,i, o)
% Notand sum(POLY1, POLY2, POLYS), semnificatia argumentelor este:
% POLY1 - primul polinom de insumat
% POLI2 - al doilea polinom de insumat
% POLYS - polinomul-suma, in forma canonica
% Insumarea se face "injectand" succesiv monoamele ambelor polinoame in acelasi
% arbore binar (vezi predicatul proper), apoi aducand rezultatul la forma canonica
sum(POLY1, POLY2, POLYS) :-
maketree(POLY1, TREE),
maketree(POLY2, TREE),
inorder(TREE, POLYSUM),
simpl(POLYSUM, POLYS),
!.
/***************************************************************/
% Calculeaza produsul a doua polinoame date in forma canonica
% (i, i, o)
% Rezultatul este intors tot sub forma canonica
prod(P1, P2, P3) :-
prodh(P1, P2, [m(0, 0)], P),
proper(P, P3), !.
/***************************************************************/
% Predicat ajutator pentru prod(poly, poly, poly)
% (i, i, o, o)
% Notand prodh(P1, P2, SUMI, PROD), semnificatia argumentelor este:
% P1 - primul polinom
% P2 - al doilea polinom
% SUMI - suma intermediara generata de inmultirile monoamelor din P1 cu P2
% PROD - la sfarsit, contine produsul P1 * P2
% Conditia de oprire:
prodh([], _, P, P).
% Produsul celor doua polinoame este produsul primului monom al
% primului polinom cu cel de-al doilea polinom, insumat cu produsul
% restului primului polinom (ignorand deci primul monom) cu al doilea polinom
prodh([MON|T], P2, SUMI, PROD) :-
prodmon(MON, P2, P),
sum(P, SUMI, SUMI1),
!,
prodh(T, P2, SUMI1, PROD).
/***************************************************************/
% Produsul unui monom cu un polinom
% (i, i, o)
% Notand prodh(MON, POLI, PROD), semnificatia argumentelor este:
% MON - monomul de inmultit
% POLI - polinomul cu care se inmulteste
% PROD - produsul MON*POLI, n forma canonica
prodmon(A, B, C) :-
prodmonh(A, B, C1),
proper(C1, C).
/***************************************************************/
% Ajutator pentru prodmon
% (i, i, o)
% Notand prodmonh(MON, POLI, PROD), semnificatia argumentelor este:
% MON - monomul de inmultit
% POLI - polinomul cu care se inmulteste
% PROD - produsul MON*POLI
% Realizeaza parcurgerea parcurgerea polinomului POLI, monom cu monom,
% inmultind pe fiecare cu monomul MON si construind lista rezultatelor
% Conditia de oprire:
prodmonh(_, [], []).
% Parcurgerea element cu element a polinomului:
prodmonh(m(C1, E1), [m(C2, E2)|T], [m(CP, EP)|REST]) :-
CP = C1 * C2,
EP = E1 + E2,
!,
prodmonh(m(C1, E1), T, REST).
/***************************************************************/
% Deriveaza un polinom:
% (i,o)
% Notand deriv(POLY1, POLY2), semnificatia argumentelor este:
% POLY1 - polinomul de derivat
% POLY2 - polinomul derivat, in forma canonica
deriv(POLY, DERIV) :-
derivh(POLY, D),
proper(D, DERIV).
/***************************************************************/
% Ajutator pentru deriv:
% (i,o)
% Notand derivh(POLY1, POLY2), semnificatia argumentelor este:
% POLY1 - polinomul de derivat
% POLY2 - polinomul derivat, care poate sa nu fie in forma canonica
% Conditia de oprire:
derivh([], []) :- !.
% Derivata unei constante este zero:
derivh([m(C, 0)|REST], RESTDERIV) :-
!,
deriv(REST, RESTDERIV).
% Derivata unui monom C*X^G este G*C*X^(G-1)
derivh([m(C, G)|REST], [m(CD, GD)|RESTDERIV]) :-
!,
CD = G * C,
GD = G - 1,
derivh(REST, RESTDERIV).
/***************************************************************/
% Aduce un polinom la forma canonica:
% (i,o)
% Notand prodh(POLI, POLIC), semnificatia argumentelor este:
% POLI - polinomul de adus la forma canonica
% POLIC - polinomul adus la forma canonica
% Aducerea la forma canonica se realizeaza prin:
% 1. "injectarea" fiecarui monom component al polinomului intr-un
% arbore binar conform relatiei "<=" aplicata gradului monoamelor
% 2. arborele binar rezultat este parcurs in inordine pentru a
% rezulta o lista de monoame ordonate crescator dupa gradul lor
% 3. lista ordonata este simplificata dupa anumite reguli, pentru
% a rezulta forma canonica
proper(POLY1, POLY2) :-
maketree(POLY1, TREE),
inorder(TREE, POLY3),
simpl(POLY3, POLY2),
!.
/***************************************************************/
% Verifica daca doua monoame sunt in relatia "<=", privitor la gradele lor
% (i,i)
% Semnificatia parametrilor este evidenta
lesseq(m(_, E1), m(_, E2)) :- E1 <= E2.
/***************************************************************/
% Transforma un polinom intr-un arbore binar, dupa regula descrisa la predicatul inject
% (i,o)
% Notand maketree(POLY, TREE), semnificatia argumentelor este:
% POLY - polinomul de transformat
% TREE - arborele binar rezultat
% injecteaza fiecare monom al polinomului in arbore, incepand din radacina:
maketree([H|T], TREE) :-
inject(TREE, H),
!,
maketree(T, TREE).
% conditia de oprire:
maketree([], _).
/***************************************************************/
% "Injecteaza" un monom intr-un arbore, dupa algoritmul:
% (m,i) - "m" inseamna ca are rol mixt
% Notand inject(NOD, MONOM), semnificatia argumentelor este:
% NOD - nodul din arborele binar la care s-a ajuns
% MONOM - monomul care se "injecteaza"
% Daca nodul este o variabila libera, "aseaza" monomul aici:
inject(X, M) :-
free(X), !, X = t(M, _, _).
% Daca gradul monomului este <= gradul nodului, propaga monomul pe ramura stanga:
inject(X, M) :-
X = t(MV, L, _),
lesseq(M, MV),
!,
inject(L, M).
% Altfel inseamna ca gradul monomului este mai mare, propaga monomul pe ramura dreapta:
inject(t(_, _, R), M) :-
!,
inject(R, M).
/***************************************************************/
% Traverseaza un arbore binar in inordine:
% (i,o)
% Notand inorder(TREE, POLY), semnificatia argumentelor este:
% TREE - arborele de traversat
% POLY - polinomul rezultat
inorder(TREE, POLY) :-
inorderh(TREE, POLY, []).
/***************************************************************/
% Ajutator pentru inorder:
% (i,o)
% Notand inorderh(TREE, POLY1, POLY2), semnificatia argumentelor este:
% TREE - nodul la care s-a ajuns
% POLY1 - elementul de lista al polinomului la care s-a ajuns acum
% POLY2 - elementul de lista liber, care va fi ocupat la urmatoarea prelucrare
% Am ajuns la o frunza-variabila libera, opreste parcurgerea:
inorderh(X, _, _) :-
free(X), !.
% Am ajuns la un nod care este un monom, insa subarborii sai sunt amandoi variabile libere.
% Ataseaza monomul din nod la lista:
inorderh(t(V, L, R), [V|X2], X2) :-
free(L), free(R),
!.
% Am ajuns la un nod-monom, cu subarborele stang var. libera, si dreptul nu.
% Ataseaza monomul la lista si continua cu parcurgerea in inordine pe subarborele drept:
inorderh(t(V, L, R), [V|X1], X2) :-
free(L), !,
inorderh(R, X1, X2).
% Am ajuns la un nod-monom, cu subarborele drept var. libera, si stangul nu.
% Parcurge in inordine subarborele stang, apoi ataseaza la lista monomul:
inorderh(t(V, L, R), X1, X2) :-
free(R), !,
inorderh(L, X1, X3),
X3 = [V|X2].
% Am ajuns la un nod-monom, care are si sub-arbore stang, si sub-arbore drept.
% Parcurge in inordine subarborele stang, ataseaza nodul, apoi parcurge sub-arborele drept:
inorderh(t(V, L, R), X1, X3) :-
!,
inorderh(L, X1, [V|X2]),
!,
inorderh(R, X2, X3).
/***************************************************************/
% Simplifica un polinom cu monoamele in ordine crescatoare a gradelor:
% (i,o)
% Notand simpl(POLY1, POLY2), semnificatia argumentelor este:
% POLY1 - polinomul cu monoamele sortate crescator ca grade
% POLY2 - polinomul in forma canonica
% Pentru cazul in care rezultatul unor operatii este un polinom vid, il egaleaza cu monomul 0*X^0:
simpl([], [m(0,0)]) :- !.
% Daca rezultatul simplificarilor este o lista vida (de ex, polinomul a fost X - X), trece la urmatoarea clauza,
% altfel se opreste aici:
simpl(POL, POLS) :-
simplh(POL, POLS), not(POLS = []), !.
% Ajunge aici doar daca rezultatul simplificarilor este un polinom-lista vida,
% si il transforma in polinomul 0*X^0
simpl(_, [m(0,0)]).
/***************************************************************/
% Ajutator pentru simpl:
% (i,o)
% Semnificatia argumentelor este aceeasi cu ale lui simpl
% Lista vida nu mai poate fi simplificata:
simplh([], []).
% Sterge un monom cu coeficient aproximativ egal cu zero:
simplh([m(Z, _)|REST], X) :-
zero(Z),
!,
simpl(REST, X).
% Aduna doua monoame de acelasi grad si continua evaluarea din acelasi punct
% (pot fi mai multe monoame de acelasi grad):
simplh([m(A, G), m(B, G)|REST], X) :-
C = A + B,
!,
simplh([m(C, G)|REST], X).
% Nu mai exista o regula de simplificare, adauga monomul la rezultat:
simplh([MON|REST], [MON|RESTS]) :-
!,
simplh(REST, RESTS).
/***************************************************************/
% PREDICATE PENTRU SCRIEREA UNUI POLINOM INTR-O FORMA LIZIBILA
% S-au folosit optimizari de genul:
% 3*X^2 se scrie 3X^2
% 1*X^N se scrie X^N
% 2*X^0 se scrie 2
writepoly([m(C, 0)]) :- !, write(C).
writepoly([m(Z, 0)|REST]) :- zero(Z), !, writepoly(REST).
writepoly([MONOM|REST]) :- writemon(MONOM), !, writepoly(REST).
writepoly([]).
writemon(m(0, _)) :- write("0"),!.
writemon(m(C, 0)) :- C<0, !, write(C).
writemon(m(C, 0)) :- C>=0, !, write("+", C).
writemon(m(UNU, 1)) :- unu(UNU), !, write("+X").
writemon(m(UNU, N)) :- unu(UNU), !,write("+X^", N).
writemon(m(MINUSUNU, 1)) :- X = - MINUSUNU, unu(X), !, write("-X").
writemon(m(MINUSUNU, N)) :- X = - MINUSUNU, unu(X), !, write("-X^", N).
writemon(m(C, 1)) :- C >0, !, write("+", C, "X").
writemon(m(C, 1)) :- C < 0, !, write(C, "X").
writemon(m(C, N)) :- C > 0, !, write("+", C, "X^", N).
writemon(m(C, N)) :- C < 0, !, write(C, "X^", N).
/***************************************************************/
% Predicat de test.
% - Cere doua polinoame
% - Le aduce la forma canonica
% - Le afiseaza
% - Calculeaza suma lor
% - Afiseaza si suma
testsuma :-
write("Test suma a 2 polinoame, P+Q."),
nl,
write("Introduceti polinomul P:"), nl,
cinpoly(P), proper(P, PP),
nl,
write("Introduceti polinomul Q:"), nl,
cinpoly(Q), proper(Q, PQ),
nl,
write("============="),
nl,
write("P(X)="), writepoly(PP), nl,
write("Q(X)="), writepoly(PQ), nl,
sum(PP, PQ, SUM),
write("P(X)+Q(X)="), writepoly(SUM), nl.
/***************************************************************/
% Predicat de test.
% - Cere doua polinoame
% - Le aduce la forma canonica
% - Le afiseaza
% - Calculeaza produsul lor
% - Afiseaza si produsul
testprod :-
write("Test produsul a 2 polinoame, P*Q."),
nl,
write("Introduceti polinomul P:"), nl,
cinpoly(P), proper(P, PP),
nl,
write("Introduceti polinomul Q:"), nl,
cinpoly(Q), proper(Q, PQ),
nl,
write("============="),
nl,
write("P(X)="), writepoly(PP), nl,
write("Q(X)="), writepoly(PQ), nl,
prod(PP, PQ, PROD),
write("P(X)*Q(X)="), writepoly(PROD), nl.
/***************************************************************/
% Predicat de test.
% - Cere un polinom
% - Il aduce la forma canonica
% - Il afiseaza
% - Calculeaza derivata sa
% - Afiseaza si derivata
testderiv :-
write("Test derivatei unui polinom P."),
nl,
write("Introduceti polinomul P:"), nl,
cinpoly(P), proper(P, PP),
nl,
write("============="),
nl,
write("P(X)="), writepoly(PP), nl,
deriv(PP, PDERIV),
write("P'(X)="), writepoly(PDERIV), nl.
/*******************************************************
* DE AICI URMEAZA SURSA INITIALA A EXERCITIULUI PROPUS:
* Nicolaie Popescu-Bodorin, 2009.
* http://fmi.spiruharet.ro/bodorin/
*******************************************************/
PREDICATES
cinmon(monom) % read a monomial in console
cinpoly(poly) % read a polynomial in console
cinpoly(integer,poly,poly)
continue(poly,integer) % end or continue to read with cinpoly
liststopoly(lr,li,poly) % from two lists to a polynomial and conversely
pow(r,i,r) % natural power of a real number
pow(r,i,r,r)
polyfeval(poly,r,r) % evaluate a polynomial function in a given point
polyfeval(poly,r,r,r)
CLAUSES
% Read a monomial in console:
cinmon(M):-
write("Give a coefficient:"),readreal(C),nl,
write("Give an exponent:"),readint(P),nl,
M = [m(C,P)],
write("Your monomial is: "),nl,
%write(C),write("X"),write(P),nl,
write(M),nl.
% Read a polynomial in console (symbolic computation):
cinpoly(P):-write("Give a polynomial."),nl,
cinpoly(1,[],P).
cinpoly(0,Ri,Ri):-
write("Your polynomial is: "),nl,
write(Ri),nl.
cinpoly(1,Ri,R):-
cinmon([M]),
continue([M|Ri],Fa),
cinpoly(Fa,[M|Ri],R).
% End a reading (cinpoly):
continue(Ri,B):-
write("Until now, your polynomial is: "),nl,
write(Ri),nl,
write("Do you want to continue?[1/0]"),readint(B).
% Convert two lists (coefficients and exponents)
% into a polynomial and conversely (symbolic computation):
% usage: (i,i,o),(o,o,i) - (lr,li,poly)
liststopoly([],[],[]).
liststopoly([HC|TC],[HE|TE],[m(HC,HE)|T]):-
liststopoly(TC, TE, T).
% Evaluate a polynomial function at a given point
% (converting a symbolic expression to a number, giving a numerical meaning to a symbolic expression):
% usage: (i,i,o) - (poly,real,real)
polyfeval([m(C,E)],X,R):-pow(X,E,P),R=C*P,!.
polyfeval(P,X,R):-polyfeval(P,X,0,R).
polyfeval([],_,Ri,Ri).
polyfeval([H|T],X,Ri,R):-
polyfeval([H],X,RH),
Ria=Ri+RH,
polyfeval(T,X,Ria,R).
% natural power of a real number:
% usage: (i,i,o) - (real,integer,real)
pow(_,0,1):-!.
pow(X,N,R):-pow(X,N,1,R).
pow(_,0,Ri,Ri).
pow(X,N,Ri,R):-N>0,Na=N-1,Ria=Ri*X,pow(X,Na,Ria,R).
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% Homework assignment:
% Implement sum, product, and derivation of polynomials.
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
]]>