Artifact 846288d7502b6314e2d0d853484b816f54629cc888facdd1eec9c4f17e0fad4e:
- Executable file
r36/src/changevr.red
— part of check-in
[f2fda60abd]
at
2011-09-02 18:13:33
on branch master
— Some historical releases purely for archival purposes
git-svn-id: https://svn.code.sf.net/p/reduce-algebra/code/trunk/historical@1375 2bfe0521-f11c-4a00-b80e-6202646ff360 (user: arthurcnorman@users.sourceforge.net, size: 9496) [annotate] [blame] [check-ins using] [more...]
module changevr; % Facility to perform CHANGE of independent VARs. %*********************************************************************; % ------------------------------- ; % C H A N G E V A R ; % ------------------------------- ; % ; % A REDUCE facility to perform CHANGE of independent VARiable(s) ; % in differential equation (or a set of them). ; % ; % ; % Author : Gokturk Ucoluk ; % Date : Oct. 1989 ; % Place : Middle East Tech. Univ., Physics Dept., Turkey. ; % Email : A07917 @ TRMETU.BITNET or A06794 @ TRMETU.BITNET ; % ; % Version: 1.00 ; % ; % ( *** Requires: REDUCE 3.0 or greater *** ) ; % ( *** Requires: Matrix package to be present *** ) ; % ; % There exists a document written in LaTeX that explains the ; % package in more detail. ; % ; % Keywords : differential equation, change of variable, Jacobian ; % ; %*********************************************************************; create!-package('(changevr),'(contrib misc)); load!-package 'matrix; fluid '(wtl!*); global '(!*derexp !*dispjacobian); switch derexp, dispjacobian ; % on derexp : Smart chain ruling % on dispjacobian : Displays inverse % Jacobian. symbolic procedure simpchangevar v; begin scalar u,f,j,dvar,flg; dvar := if pairp car v then cdar v else car v.nil; v := cdr v; % Dvar: list of depend. var u := if pairp car v then cdar v else car v.nil; v := cdr v; % u: list of new variables if eqcar(car v,'list) then v := append(cdar v,cdr v); while cdr v do << if caar v neq 'equal then rederr "improper new variable declaration"; f := cdar v . f; % f: list of entries (oldvar func(newvrbs)) v := cdr v >>; % i i v := reval car v; % v: holds now last expression (maybe a list) if length u < length f then rederr "Too few new variables" else if length u > length f then rederr "Too few old variables"; % Now we form the Jacobian matrix ; j := for each entry in f collect for each newvrb in u collect reval list('df,cadr entry, newvrb); j := cdr aeval list('quotient,1,'mat.j); % j: holds inverse Jacobian. % We have to define the dependencies of old variables to new % variables. for each new in u do for each old in f do depend1(new, car old, t); % Below everything is perplexed : % The aim is to introduce LET DF(new ,old ) = jacobian % row col row,col % With the pairing trick below we do it in one step. % new : car row, old : caar col, jacobian : cdr col % row col row,col % for each row in pair(u,j) do for each col in pair(f,cdr row) do << let2(list('df,car row,caar col), sqchk cdr col, nil, t); if !*dispjacobian and !*msg then mathprint list('equal,list('df,car row,caar col), sqchk cdr col) >>; flg := !*derexp; !*derexp := t; v := changearg(dvar,u,v); for each entry in f do v := subcare(car entry, cadr entry, v); % now here comes the striking point ... we evaluate the last % argument. v := simp!* v; % Now clean up the mess of LET; for each new in u do for each old in f do << let2(list('df,new,car old), nil, nil, nil); let2(list('df,new,car old), nil, t, nil) >>; !*derexp := flg; return v; end; put('changevar,'simpfn,'simpchangevar); symbolic procedure changearg(f,u,x); if atom x then x else if car x memq f then car x . u else changearg(f,u,car x) . changearg(f,u,cdr x); symbolic procedure subcare(x,y,z); % shall be used after changearg ; if null z then nil else if x = z then y else if atom z or get(car z,'subfunc) then z else (subcare(x,y,car z) . subcare(x,y,cdr z)); % Updated version of DIFFP.. chain rule handling is smarter. ; % Example: If F is an operator and R has a implicit dependency on X, % declared by a preceding DEPEND R,X .. then the former version % of DIFFP, provided in REDUCE 3.3, was such that an algebraic % evaluation of DF(F(R),X) will evaluate to itself, that % means no change will happen. With the below given update this % is improved. If the new provided flag DEREXP is OFF then % the differentiation functions exactly like it was before, % but if DEREXP is ON then the chain rule is taken further to % yield the result: DF(F(R),R)*DF(R,X) . remflag('(diffp),'lose); % Since we want to reload it. symbolic procedure diffp(u,v); %U is a standard power, V a kernel. % Value is the standard quotient derivative of U wrt V. begin scalar n,w,x,y,z,key; integer m; n := cdr u; %integer power; u := car u; %main variable; if u eq v and (w := 1 ./ 1) then go to e else if atom u then go to f %else if (x := assoc(u,dsubl!*)) and (x := atsoc(v,cdr x)) % and (w := cdr x) then go to e %deriv known; %DSUBL!* not used for now; else if (not atom car u and (w:= difff(u,v))) or (car u eq '!*sq and (w:= diffsq(cadr u,v))) then go to c %extended kernel found; else if x := get(car u,'dfform) then return apply3(x,u,v,n) else if x:= get(car u,dfn_prop u) then nil else if car u eq 'plus and (w:=diffsq(simp u,v)) then go to c else go to h; %unknown derivative; y := x; z := cdr u; a: w := diffsq(simp car z,v) . w; if caar w and null car y then go to h; %unknown deriv; y := cdr y; z := cdr z; if z and y then go to a else if z or y then go to h; %arguments do not match; y := reverse w; z := cdr u; w := nil ./ 1; b: %computation of kernel derivative; if caar y then w := addsq(multsq(car y,simp subla(pair(caar x,z), cdar x)), w); x := cdr x; y := cdr y; if y then go to b; c: %save calculated deriv in case it is used again; %if x := atsoc(u,dsubl!*) then go to d %else x := u . nil; %dsubl!* := x . dsubl!*; % d: rplacd(x,xadd(v . w,cdr x,t)); e: %allowance for power; %first check to see if kernel has weight; if (x := atsoc(u,wtl!*)) then w := multpq('k!* .** (-cdr x),w); m := n-1; % Evaluation is far more efficient if results are rationalized. return rationalizesq if n=1 then w else if flagp(dmode!*,'convert) and null(n := int!-equiv!-chk apply1(get(dmode!*,'i2d),n)) then nil ./ 1 else multsq(!*t2q((u .** m) .* n),w); f: % Check for possible unused substitution rule. if not depends(u,v) and (not (x:= atsoc(u,powlis!*)) or not smember(v,simp cadddr x)) then return nil ./ 1; w := list('df,u,v); go to j; h: %final check for possible kernel deriv; y := nil; if car u eq 'df then key:=t; w := if key then 'df . cadr u . derad(v,cddr u) else list('df,u,v); y := cddr u; w := if (x := opmtch w) then simp x else if not depends(cadr w,caddr w) then nil ./ 1 else if !*derexp then begin if atom cadr w then return mksq(w,1); w := nil ./ 1; for each m in cdr(if key then cadr u else u) do w := addsq(multsq( if (x := opmtch (z := 'df . if key then (cadr u.derad(m,y)) else list(u,m) )) then simp x else mksq(z,1), diffsq(simp m,v)), w); return w end else mksq(w,1); go to e; j: w := if x := opmtch w then simp x else mksq(w,1); go to e end; endmodule; end;