(*^

::[paletteColors = 128; 
	fontset = title, "New York", 24, L3, center, bold, nohscroll;
	fontset = subtitle, "New York", 18, L2, center, bold, nohscroll;
	fontset = subsubtitle, "New York", 14, L2, center, bold, nohscroll;
	fontset = section, "New York", 14, L2, bold, nohscroll, grayBox;
	fontset = subsection, "New York", 12, L2, bold, nohscroll, blackBox;
	fontset = subsubsection, "New York", 10, L2, bold, nohscroll, whiteBox;
	fontset = text, "New York", 12, L2, nohscroll;
	fontset = smalltext, "New York", 10, L2, nohscroll;
	fontset = input, "Courier", 12, L2, bold, nowordwrap;
	fontset = output, "Courier", 12, L2, nowordwrap;
	fontset = message, "Courier", 12, L2, R65535, nowordwrap;
	fontset = print, "Courier", 12, L2, nowordwrap;
	fontset = info, "Courier", 12, L2, nowordwrap;
	fontset = postscript, "Courier", 12, L2, nowordwrap;
	fontset = name, "Geneva", 10, L2, italic, B65535, nowordwrap, nohscroll;
	fontset = header, "Times", 10, L2;
	fontset = footer, "Times", 12, L2, center;
	fontset = help, "Geneva", 10, L2, nohscroll;
	fontset = clipboard, "New York", 12, L2;
	fontset = completions, "New York", 12, L2, nowordwrap;
	fontset = network, "Courier", 10, L2, nowordwrap;
	fontset = graphlabel, "Courier", 12, L2, nowordwrap;
	fontset = special1, "New York", 12, L2, nowordwrap;
	fontset = special2, "New York", 12, L2, center, nowordwrap;
	fontset = special3, "New York", 12, L2, right, nowordwrap;
	fontset = special4, "New York", 12, L2, nowordwrap;
	fontset = special5, "New York", 12, L2, nowordwrap;]
:[font = input; initialization; ]
*)
(* Function and variable differential substitutions.
   Should be working but is not guaranteed to be bug-free!
   PD software 4.4.1990 by
           Pekka Janhunen
           Finnish Meteorological Institute
           Geophysics Dept.
           pjanhune@finsun.csc.fi  *)

BeginPackage["DSubst`"]

FuncChange::usage = "FuncChange[expr, y, u, f]
  performs a change of function y = f(u) in expr.
  Example: FuncChange[ y[x]+y'[x], y,u,Exp[u] ] returns 
  Exp[u[x]] + Exp[u[x]] u'[x]. FuncChange[ y1'[x] + y2'[x],
  {y1,y2}, {u1,u2}, {Exp[u1], 2 u2} ] returns Exp[u1[x]] u1'[x] + 2 u2'[x]."
  
FuncChange::notes = "FuncChange calls FCollect after
  finishing its work. The options for FCollect can be
  given to FuncChange, too."

VarChange::usage = "VarChange[expr, y,x,t,f] performs a change of variables
  x = f(t) in expr when y is dependent function.
  Example: VarChange[y'[x] + y[x], y, x, t, 2t] returns
  y[t] + y'[t]/2."
  
FCollect::usage = "FCollect[F,y,{options}] collects F according to y."

SimplifyEqn::usage = "SimplifyEqn[F,x] removes factors from F which are
  != 0 for all x."

Begin["`Private`"]

(* Map3 = LISP's MAPCAR-routine for 3 arguments: *)

Map3[fn_, args1_List, args2_List, args3_List]:=
  Prepend[Map3[fn,Rest[args1],Rest[args2],Rest[args3]],
          Apply[fn,List[First[args1],First[args2],First[args3]]]];

Map3[fn_,args1_List, args2_List, args3_List]:=
  List[Apply[fn, List[First[args1], First[args2], First[args3]]]] /; Length[args1] == 1

Map3[fn_, {}, {}, {}]:= {}



FuncChange[F_, y_List, u_List, f_List, opts___]:=
  Block[ {rules = {}},
         Map3[ AppendTo[rules, #1[x__] :> (#3/.#2->#2[x])]&,
               y,u,f ];
         Map3[ AppendTo[rules,
                        Derivative[n__][#1][x__] :>
                          FuncChangeD[(#3/.#2->#2[x]), {x}, {n} ]]&,
               y,u,f];
         FCollect[F/.rules, u, opts] ]  /; Length[y] == Length[u] == Length[f]

FuncChange[F_, y_Symbol, u_Symbol, f_, opts___]:= FuncChange[F,{y},{u},{f}, opts]

FuncChangeD[f_,x_List,n_List]:=
  Apply[D,Prepend[Transpose[{x,n}],f]] /; Length[x] == Length[n]



JacobianMatrix[x_List,t_List]:=
  Block[ {i,j},
         Table[ D[x[[i]],t[[j]]], {i,Length[x]}, {j,Length[t]}]] /;
  Length[x] == Length[t]

DOper[n_List]:=
  Block[ {i, result = {} },
         Do[ Do[AppendTo[result, D1[i]], {n[[i]]}],{i,Length[n]}];
         Apply[HeadCompose,result] ]

VarChange[F_, yy_List, x_List, t_List, f_List]:=
  (Block[{InverseJacobian = Inverse[JacobianMatrix[f,t]]},
         F /.
         {Apply[Derivative[n__][y_?(MemberQ[yy,#]&)],x] :>
             (DOper[{n}][Apply[y,t]] /. D1[i_] :> InverseJacobian[[i]].Array[d,Length[x]]
                        /. d[i_] :> (D[#,t[[i]]]&) ),
          Apply[y_?(MemberQ[yy,#]&),x] :> Apply[y,t]}
                         ] //.
   Times[anything1_,Function[anything2_]]:>Function[anything1 anything2]
   /. x1_?( (AtomQ[#] && MemberQ[x,#] && Length[DSubst`pos = Position[x,#]] > 0) &) :> 
      f[[DSubst`pos[[1,1]]]])

VarChange[F_, y_Symbol, x_Symbol, t_Symbol, f_]:= VarChange[F,{y},{x},{t},{f}]


Options[FCollect]:= {Prolog -> Simplify}

FCollect[F_, y_List, opts___]:= 
  Block[ {Dys = {}},
         F /. {y1_Symbol?(MemberQ[y,#]&)[x__] :>
                 (AppendTo[Dys, y1[x]];Null),
               Derivative[n__][y1_Symbol?(MemberQ[y,#]&)][x__] :>
                 (AppendTo[Dys, Derivative[n][y1][x]];Null) };
         Dys = Union[Dys];
         Collect[F,Dys] /.
           anything_ yy_?(MemberQ[Dys,#]&) :>
             (Prolog/.{opts}/.Options[FCollect])[anything] yy]

FCollect[F_,y_Symbol, opts___]:= FCollect[F,{y},opts]

SimplifyEqn[eqn_,x_]:=
  Apply[Times,Cases[ Apply[List,FactorTerms[eqn,x]], !TrueQ[(# == 0) == False] ]]

End[]

Print[{FuncChange, VarChange, FCollect, SimplifyEqn}]

EndPackage[]

Null





(* The following work but are unnecessary:

FuncChange1[F_,y_,u_,f_]:=
  (F /.
   {Derivative[n_][y][x_] :> D[f/.u->u[x],{x,n}],
    y[x_] :> f/.u->u[x]} )

VarChange1[F_,x_,t_,f_]:=
  (F /.
   {Derivative[n_Integer][y_][x] :> Nest[(D[#,t]/D[f,t])&, y[t],n],
    y_[x] :>  y[t] } )

FindArgs[f_,y_]:=
  Block[ {pos = Union[Position[f,y[x__]],Position[f,Derivative[n__][y][x__]]],
          result={}},
         MapAt[((result = Union[result,#/.Apply[List,#]])&),f,pos];
         result ]*)
(*
^*)