```%	diff.rb
%
%	Purpose:	(Partial) differentiation
%	Copyright:	R.A. van Engelen, Leiden University, 1996
%
%	Remarks:	Only to be used in combination with diff script!

include linear.

default strategy diff(full).
full strategy [linex, diff, opdiff, chint, fold].
chain strategy [chdiff, opdiff, chint, fold].
none strategy [nodiff].

diff; use_module([depend]).

(diff; diff(Var, X) => DYX <- is_var(Var), ( (to_atom(X, D), D has dependency(Zs) ->    map(diff(Var, Zs), DVZs),    map(diff(Zs, X), DZsX),    map(DVZs * DZsX, DVZX),    DYX = (+ | DVZX)   ) ; DYX = df(Var, X) )).
diff; diff(Y1 + Y2, X) => diff(Y1, X) + diff(Y2, X).
diff; diff(- Y, X) => - diff(Y, X).
diff; diff(undefined, _G596) => {undefined}.
diff; diff(1 / Y, X) => diff(Y, X) * Y ^ {- {2}}.
diff; diff((* | [Y | Ys]), X) => ( (* | [diff(Y, X) | Ys])                                 + Y * diff((* | Ys), X)                                 ).
diff; diff(Y ^ N, X) => ( N                        * diff(Y, X)                        * Y ^ (N - 1)                        ) <- is_number(N).
diff; diff(Y ^ Z, X) => diff(exp(log(Y) * Z), X).
diff; diff(prod(Y, D), X) => sum(diff(Y, X) / Y * prod(Y, D), D).
diff; diff(X, X) => {1}.

%  d(if U then V else W ) = if U then d(V) else d(W).
diff; diff(A if B, X) => diff(A, X) if B.
diff; diff(A \\ B, X) => diff(A, X) \\ diff(B, X).

%  stencil activities
diff; diff(protect(A), X) => protect(diff(A, X)).
diff; diff(protect(X, Y), Z) => protect(diff(X, Z), Y).

diff; diff(varidx(X, _G596), X) => {1}.
diff; diff(varidx(Var, Idx1), varidx(Var, Idx2)) => if((and | EqIdx), {1}, {0}) <- map(Idx1 == Idx2, EqIdx).
(diff; diff(varidx(Var, Idx), X) => + | DUIdxX <- X = {Y}, Y instance_of index_var, Var has point(Pnt), findall(diff(J, X) * df({varidx({Var}, IdxI)}, I), ( nth1(N, Pnt, I)                                                    , nth1(N, Idx, J)                                                    , replace_nth1(N, Idx, I, IdxI)                                                    ), DUIdxX)).
% independent variables behave like constants, so linearise to 0
diff; diff(I, nil) => 0 <- to_atom(I, AI), AI instance_of index_var.
(diff; diff(L, nil) => 0 <- to_atom(L, AL), AL instance_of linearized_var).
(diff; diff(L, nil) => 0 <- to_atom(L, AL), AL instance_of parameter_var).

% return the linearised variable, linearise the discrete value
(diff; diff(U, nil) => DU <- to_atom(U, AU), AU instance_of variable, !, ( AU has linearized_var(DU) ; ( (concat(AU, '_TL', DU), \+ (DU isa object) ->      assume DU isa linearized_var     )   ; linearized_var do create(DU)   ),   AU do get_properties(Ps),   assume DU has Ps,   ignore(remove DU has discrete_value(_G653)),   ignore(remove DU has fortran_name(_G662)),   assume AU has linearized_var(DU),   assume DU has parent_var(AU) ), ( AU instance_of fundamental_var ; DU has discrete_value(_G709) ; (AU has discrete_value(DAU) ->    normalize(diff(DAU, nil), NDDAU),    denormalize(NDDAU, DDAU),    assume DU has discrete_value(DDAU)   ) )).
(diff; diff(S, X) => OpD <- (to_op(S, Op) ->  Op instance_of stencil_op,  Op =.. [F, A | As],  OpD =.. [F, diff(A, X) | As] )).
% needed: h_x is inert but that should be clear??!??
% worse: kappa is inert but that is not accepted: assigned!
% chain rules
(diff; diff({varidx({Var}, Idx)}, X) => DAX <- to_atom(Var, A), to_atom(X, D), A has dependency(Zs), ( memberchk(D, Zs) -> DAX = df(varidx({Var}, Idx), X) ; findall(diff(varidx({Var}, Idx), Z) * diff(Z, X), member(Z, Zs), DAZX), DAX = (+ | DAZX) ), !).
(diff; diff(Atom, X) => DAX <- to_atom(Atom, A), to_atom(X, D), ( (A has dependency(Zs) ->    ( memberchk(D, Zs) -> DAX = df(A, X)    ; findall(diff(A, Z) * diff(Z, X), member(Z, Zs), DAZX),      DAX = (+ | DAZX)    )   ) ; D has dependency(Zs), memberchk(A, Zs) -> DAX = df(A, X) ; DAX = {0} ), !).

diff; diff(Const, _G596) => {0} <- is_constant(Const).
(diff; diff(Op, X) *=> + | DFAsX <- to_op(Op, FAs), !, \+ (FAs do commute_over(diff(_G610, X))), \+ (diff(_G620, X) do commute_over(FAs)), findall(diff(A, X) * df(FAs, N), arg(N, FAs, A), DFAsX)).

%	Symbolic differentiation by chain-rule only

chdiff; diff(undefined, _G596) => undefined.
chdiff; diff(X, X) => {1}.
(chdiff; diff(X, Y) => DXY <- to_atom(X, A), \+ (A has dependency(_G607)), to_atom(Y, B), ( (B has dependency(Zs) ->    (memberchk(A, Zs) -> DXY = df(X, Y); DXY = {0})   ) ; DXY = {0} )).
(chdiff; diff(FAs, X) => DYX <- to_atom(X, D), ( (depend : infer(FAs, [dep(Zs)]) ->    ( memberchk(D, Zs) -> DYX = df(FAs, X)    ; findall(diff(FAs, Z) * diff(Z, X), member(Z, Zs), DYZX),      DYX = (+ | DYZX)    )   ) ; (D has dependency(Zs) ->    map(diff(FAs, Zs), DFZs),    map(diff(Zs, X), DZsX),    map(DFZs * DZsX, DYZX),    DYX = (+ | DYZX)   ) )).
chdiff; diff(E, X) => df(E, X).

chint; integrate(undefined, _G596) => undefined.
chint; integrate(X, X = L .. U) => 0 max ( 1 / 2                                         * (U ^ 2 - L ^ 2)                                         ).
(chint; integrate(X, Y = L .. U) => IXY <- to_atom(X, A), \+ (A has dependency(_G613)), to_atom(Y, B), B has dependency(Zs), ( memberchk(A, Zs) -> IXY = int(X, Y = L .. U) ; IXY = X * (U - L) )).
(chint; integrate(FAs, X = L .. U) => IYX <- depend : infer(FAs, [dep(Zs)]), to_atom(X, D), ( memberchk(D, Zs) -> IYX = int(FAs, X = L .. U) ; findall(integrate(FAs * diff(X, Z), Z = Z \$ L .. Z \$ U), member(Z, Zs), IYZX), IYX = (+ | IYZX) ), !).
chint; integrate(X, Y) => int(X, Y).

%	No symbolic differentiation/integration

nodiff; diff(E, X) => df(E, X).
nodiff; integrate(E, X) => int(E, X).

%	Functions that are not differentiable on whole domain:

opdiff; df({abs(X)}, {1}) => if(X <> {0}, sign(X), {undefined}).
opdiff; df({ceil(X)}, {1}) => if(X <> ceil(X), {0}, {undefined}).
opdiff; df({floor(X)}, {1}) => if(X <> floor(X), {0}, {undefined}).
opdiff; df({frac(X)}, {1}) => if(X <> trunc(X), {1}, {undefined}).
opdiff; df({round(X)}, {1}) => if(X <> round(X), {0}, {undefined}).
opdiff; df({trunc(X)}, {1}) => if(X <> trunc(X), {0}, {undefined}).

%	Functions that are differentiable:

opdiff; df(acos(X), 1) => - ({1} - X ^ {2}) ^ {- {{1} / {2}}}.
opdiff; df(acosh(X), 1) => (X ^ {2} - {1}) ^ {- {{1} / {2}}}.
opdiff; df(asin(X), 1) => ({1} - X ^ {2}) ^ {- {{1} / {2}}}.
opdiff; df(asinh(X), 1) => ({1} + X ^ {2}) ^ {- {{1} / {2}}}.
opdiff; df(atan(X), 1) => {1} / ({1} + X ^ {2}).
opdiff; df(atanh(X), 1) => {1} / ({1} - X ^ {2}).
opdiff; df(cos(X), 1) => - sin(X).
opdiff; df(cosh(X), 1) => - sinh(X).
opdiff; df(exp(X), 1) => exp(X).
opdiff; df(log(X), 1) => 1 / X.
opdiff; df(sign(_G595), 1) => 0.0.
opdiff; df(sin(X), 1) => cos(X).
opdiff; df(sinh(X), 1) => cosh(X).
opdiff; df(sqrt(X), 1) => {{1} / {2}} * X ^ {- {{1} / {2}}}.
opdiff; df(tan(X), 1) => {1} + tan(X) ^ {2}.
opdiff; df(tanh(X), 1) => {1} - tanh(X) ^ {2}.
(opdiff; df(Op, N) => DOp <- to_integer(N, I), to_op(Op, FAs) -> FAs has derivative(I, DOp)).

```