%********************************************************************
module decoupling$
%********************************************************************
%  Routines for decoupling de's
%  Author: Andreas Brand 1995, updates and extensions by Thomas Wolf
%  Comments on some routines by Arrigo.
%
%  $Id: crdec.red,v 1.10 1998/06/25 18:19:36 tw Exp tw $
%

symbolic procedure which_deriv(p,q)$
%  yields a list of variables and orders 
%  such that one gets at least q by differentiating p w.r.t. the vars 
%  p,q: lists of variables and orders
begin scalar l,n,a$
while q do
   if (a:=member(car q,p)) then
      <<q:=cdr q$
      if q and numberp(car q) then
      	 <<n:= car q$
         q:=cdr q>>
      else n:=1$
      n:=n-(if pairp cdr a and numberp cadr a then cadr a else 1)$
      if n>0 then
      	 <<l:=cons(car a,l)$
         if n>1 then l:=cons(n,l)>> >>
   else 
      <<l:=cons(car q,l)$
      q:=cdr q$
      if q and numberp(car q) then 
	 <<l:=cons(car q,l)$
	 q:=cdr q>> >>$
return append(reverse l,q)$
end$

symbolic procedure dec_ld_info(p,q,simpp,simpq,f,vl,rl)$
%  gets leading derivatives of f in p and q wrt. vars order vl
%  and the lists of variables and orders for differentiation
begin scalar s,l,l1,l2,vl1,vl2,d,ld1,ld2,wd1,wd2$
 %
 % if (p has more variables than q) or 
 %    (f is not leading function of p)
 % => simpp = t => p must be simplified with (deriv.s of) q 
 % if (q has more variables than p) or
 %    (f is not leading function of q)
 % => simpq = t => q must be simplified with (deriv.s of) p 
 %
 % vl1 holds the list of _ordered_ variables of p
 % vl2 holds the list of _ordered_ variables of q
 %
 % list all powers of derivatives of f in p as l1 and in q as l2
 %

 if simpp and simpq then return nil$

 vl1:=intersection(vl,get(p,'vars))$
 vl2:=intersection(vl,get(q,'vars))$

 % drop the powers from l1
 for each a in get(p,'derivs) do
  if caar a=f then l1:=cons(a,l1)$
 l1:=sort_derivs(reverse l1,list f,vl1)$
 %
 % l1 is a list of _all_ derivatives of f in p _sorted_ stored as a
 % dotted pair, e.g. ((f x 2 y) . 5) would be f_{xxy}^5, or more
 % generally ((f_1 . power) (f_2 . power) ... )
 %
 %terpri()$write "l10=",l1$ 
 l:=nil$
 for each a in l1 do
  if not member(cdar a,l) then l:=cons(cdar a,l)$
 % 
 % cdar a is the function in the pair so we are making sure that our
 % list l1 has no repetitions
 %
 l1:=reverse l$ % e.g. l1 = ((x 2 y) (x y 2) ...)
 %terpri()$write "l11=",l1$ 

 %
 % The above now applies but with q and l2 instead of p and l1
 %
 % drop the powers from l2
 for each a in get(q,'derivs) do
  if caar a=f then l2:=cons(a,l2)$
 l2:=sort_derivs(reverse l2,list f,vl2)$
 %terpri()$write "l20=",l2$ 
 l:=nil$
 for each a in l2 do
  if not member(cdar a,l) then l:=cons(cdar a,l)$
 l2:=reverse l$ % e.g. l2 = ((x 2 y) (x y 2) ...)
 %terpri()$write "l21=",l2$ 

 % At this point we have two lists, l1 and l2 resp. containing the
 % sorted list of all derivatives of the fn f in p and q and we
 % calculate the absolute degree of the highest derivative in l1

 if not simpp then << 
   %p may be differentiated and q be substituted or new equ. added 
   d:=absdeg(car l1)$
   %terpri()$write "d=",d$
   l:=l2$
   while l and (d<=absdeg car l) do
   <<s:=which_deriv(car l1,car l)$
     %terpri()$write "s1=",s$
     %
     % which_deriv(a,b) takes two lists of derivatives and returns how
     % often you need to diff. a in order to get at least the
     % derivatives in b.
     % e.g. which_deriv((x 2 y), (x y 2)) returns y
     %
     if (absdeg s + d)=absdeg car l then 
       %
       % At this point we compare the degree of the highest
       % derivative of l1 + number of diff. in order to get the
       % leading deriv. of l2 (aliased to l)
       %
       <<ld2:=car l$ l:=nil>>
     else l:=cdr l
   >> 
 >>$

 % Now, either l is nil and ld2 = leading deriv. of l2 (i.e. highest
 % deriv. of in q) [case in which leading deriv. in l2 can be obtained
 % by diff. of the leading deriv. in l1] OR ld2 is nil and l contains
 % the rest of the deriv. of l2 except the leading one [in this case
 % we _cannot_ obtain the leading deriv. in l2 by diff. the leading
 % deriv. in l1].

 if (not ld2) and (not simpq) then <<
   %
   % We cannot get to the leading deriv. in l2 by diff. of leading
   % deriv. in l1.
   %
   d:=absdeg(car l2)$
   %
   % We try the opposite way, we try to diff. something in l2 to get
   % into l1.
   %
   l:=l1$
   while l and (d<=absdeg car l) do
   <<s:=which_deriv(car l2,car l)$
     %terpri()$write "s2=",s$ 
     if (absdeg s + d)=absdeg car l then <<ld1:=car l$ l:=nil>>
                                    else l:=cdr l
   >> 
 >>$

 % We now have either ld2 non-nil, i.e. we can get to leading derv. in
 % l2 by differentiation of terms in l1 OR we have ld1 non-nil in
 % which case we have the opposite situation. If neither are non-nil
 % then we have to cross-differentiate to get the ld to match.
 %
 % What we return is
 %
 % ( (s ld(l1)) (nil ld(l2)) )		[ld2 non-nil] or
 % ( (nil ld(l1)) (s ld(l2)) )		[ld1 non-nil] or
 % ( (v ld(l1)) (w ld(l2)) )		[both ld1 _and_ ld2 nil]
 %
 % where v and w are the required diff. to get to ld2 and ld1 resp.
 % and s is the required diff. for the non-nil cases.
 %
 % It is to be interpreted as:
 %
 % Either	"diff. ld(l1) by s to get ld(l2)" or
 %		"diff. ld(l2) by s to get ld(l1)" or
 %		"diff. ld(l1) by wd1 and ld(l2) by wd2 to get the 
 %		ld's to match".
 %

 return
 if ld2 then list(list(s,car l1),list(nil,ld2)) else 
 if ld1 then list(list(nil,ld1),list(s,car l2)) else <<
   wd1:=which_deriv(car l1,car l2)$
   wd2:=which_deriv(car l2,car l1)$
   if (simpq and wd2) or 
      (simpp and wd1) or
      (rl and wd1 and wd2) then nil 
   else list(list(wd1,car l1),
             list(wd2,car l2))
 >>
end$

symbolic procedure dec_info(p,q,f,vl,rl)$
% yields information for a decouple reduction step
% i.e. ((info_1 info_2 deriv_to_eliminate equ_to_be_subst).num_value)
% where num_value is a measure of cost, e.g.
% (((e4 (x 2 y)) (e5 (z)) (df f x 2 y 2 z) nil) . num_value)
% Criteria:   a) the function f must depend on all vars
%             b) the function and all their derivatives must occur
%                polynomially
begin scalar a,b,l,l1,info,m,fp,fq,fpmq,fqmp,s,lenp,lenq$
  %
  % 'length is the property containing the expression length
  %
  lenp:=get(p,'length)$
  lenq:=get(q,'length)$
  if rl and ((lenp*lenq)>max_red_len) then return nil;

  a:=get(p,'vars); b:=get(q,'vars);
  l:=dec_ld_info(p,q,
                 (null get(p,'allvarfcts)) or    % star-equn.
                 (f neq car get(p,'allvarfcts)), % f not lead. fct.
                 (null get(q,'allvarfcts)) or
                 (f neq car get(q,'allvarfcts)),
                 f,vl,rl)$
  if not l then <<
    add_both_dec_with(f,p,q,rl)$
    return nil
  >>$
  %
  % l:= dec_ld_info(p,q,f,vl,rl) returns a list of lists, from which
  % a := caar l sets a to be the differentiations required to get
  % the ld(p) w.r.t. f to match that of ld(q) w.r.t. f,
  % b := caadr l sets b to be the differentiations required to get
  % the ld(q) w.r.t. f to match that of ld(q) w.r.t. f.
  %
  % l1 := cadadr l sets l1 to be the derivative in q which we
  % eliminate, similarly l is the derivative in p which we elim.
  %
  a:=caar l$             % a are the differentiations of p
  b:=caadr l$            % b are the differentiations of q
  l1:=cadadr l$
  l:=cadar l$
  % s is the equation to be substituted
  if a and not b then s:=q                 % p will be diff.
  else if b and not a then s:=p            % q will be diff.
  else if not (a or b) then                % no diff., only reduction
   if car get(p,'fcts) neq f then s:=p else  % no new leading fncts
   if car get(q,'fcts) neq f then s:=q else
   if car get(p,'derivs)=car get(q,'derivs) then
     if lenp<=lenq then s:=q
                   else s:=p                else 
%   if lower_lderiv(p,q,union(get(p,'fcts),get(q,'fcts)),vl) then s:=q
%                                                            else s:=p$
     if abs_dfrel(car get(p,'derivs),car get(q,'derivs),vl) then s:=q
                                                            else s:=p$

  fp:=get(p,'allvarfcts)$ % functions of all vars in p
  fq:=get(q,'allvarfcts)$ % functions of all vars in q

  % If a pde will be replaced by a pde with more fcts of all vars
  % then this pairing will have a lowered priority
  fqmp:=length setdiff(fq,fp);
  fpmq:=length setdiff(fp,fq);
%if nil then
  if tr_decouple then << terpri()$
    write"p=",p," q=",q," s=",s," lfp=",length fp,
         " lfq=",length fq," lfu=",length union(fp,fq),
         " fqmp=",fqmp," fpmq=",fpmq
  >>$
  m:=(1.5^absdeg(a)*lenp+1.5^absdeg(b)*lenq)*
     (length union(fp,fq))**3$ 
%if nil then
  if tr_decouple then write" m2=",m;
  if s then << 
    % the equation s will be replaced by the new one
    if fcteval(s,nil) then m:=m*10**7; 
    if (s=q) and (lenp>lenq) then m:=(m*lenp)/lenq else
    if (s=p) and (lenq>lenp) then m:=(m*lenq)/lenp;
    if (s=p) and (fqmp>0) then m:=m*10**(2*fqmp) else
    if (s=q) and (fpmq>0) then m:=m*10**(2*fpmq)
  >>   else 
  % Enlarge m because extra equation is generated (temp. idea)
  m:=m*10$
  % Non-linearity in largest derivative not taken care of.
%if nil then
  if tr_decouple then write" m3=",m;
  info:=cons(list(list(p,a),list(q,b),
                  reval cons('DF,cons(f,append(l,a))),s),m)$
  return info$
end$

symbolic procedure dec_put_info(l,rl)$
% l has form ((e4 (x 2 y)) (e5 (z)) (df f x 2 y 2 z) nil)
% puts informations for decouple reduction step
% result: ((df f x 2 y 2 z) e4 e5 nil)
if l then
begin scalar f$
  put(caar l,'dec_info,cadar l)$     % saves (x 2 y) for e4
  put(caadr l,'dec_info,cadadr l)$   % saves (z)     for e5
  if (cadar l) and (cadadr l) then << % if both eq. are diff.
    f:=caddr l;
    if pairp f then f:=cadr f;
    add_both_dec_with(f,caar l,caadr l,rl)$
  >>$
  return list(caddr l,caar l,caadr l,cadddr l)$
end$

% symbolic procedure dec_put_info(f,l)$
% % put informations for decouple reduction step
% % result: (deriv_to_eliminate pde_1 pde_2)
% if l then
% begin scalar a,b$
%    put(caar l,'dec_info,cadar l)$
%    a:=get(caar l,'dec_with)$
%    b:=assoc(f,a)$
%    a:=delete(b,a)$
%    if b then b:=cons(f,cons(caadr l,cdr b))
%         else b:=list(f,caadr l)$
%    put(caar l,'dec_with,cons(b,a))$
%    put(caadr l,'dec_info,cadadr l)$
%    a:=get(caadr l,'dec_with)$
%    b:=assoc(f,a)$
%    a:=delete(b,a)$
%    if b then b:=cons(f,cons(caar l,cdr b))
%         else b:=list(f,caar l)$
%    put(caadr l,'dec_with,cons(b,a))$
% return list(caddr l,caar l,caadr l)$
% end$

%% symbolic procedure dec_info_leq(p,q)$
%% %  relation "<=" for decouple informations
%% if p and q then 
%%    if not (cadar car p and cadadr car p) then 
%%       if not (cadar car q and cadadr car q) then (cdr p<=cdr q)
%%                                     else p
%%    else if cadar car q and cadadr car q then (cdr p<=cdr q)
%%                                 else nil
%% else p$ 

symbolic procedure dec_info_leq(p,q)$
%  relation "<=" for decouple informations
if p and q then (cdr p<=cdr q)
else if p then p
else q$

symbolic procedure dec_and_fct_select(pdes,vl,rl)$
% select 2 pdes for decoupling
% if rl then one pde must be simplified with the help of 
% another one and reduce its length
begin scalar min,f,l,l1,l2,done_pdes,car_pdes$
  while pdes do <<
    car_pdes:=car pdes;
    if flagp(car_pdes,'to_decoup) and 
       get  (car_pdes,'allvarfcts) then
    <<f:=car get(car_pdes,'allvarfcts)$
      l:=assoc(f,get(car_pdes,'dec_with))$
      if rl then l:=append(l,assoc(f,get(car_pdes,'dec_with_rl)))$
      %--- unchecked pairings 
      for each p in cdr pdes do
      if flagp(p,'to_decoup) and 
         member(f,get(p,'rational)) and
         not member(p,l) and
         not member(car_pdes,assoc(f,get(p,'dec_with))) 
      then <<l1:=dec_info(car_pdes,p,f,vl,rl)$
             if l1 then l2:=cons(l1,l2)
           >>$
      %--- checked pairings but where f is not the
      %--- leading function in car done_pdes which has not
      %--- been checked when this pairing was tested before
      for each p in done_pdes do
      if flagp(p,'to_decoup) and 
         member(f,get(p,'rational)) and
         not member(p,l) and
         not member(car_pdes,assoc(f,get(p,'dec_with))) and
         ((null get(p,'allvarfcts)) or
          (f neq car get(p,'allvarfcts)))  
      then <<l1:=dec_info(car_pdes,p,f,vl,rl)$
             if l1 then l2:=cons(l1,l2)
           >>$
    >>$
    done_pdes:=cons(car_pdes,done_pdes)$
    pdes:=cdr pdes
  >>$

  %--- l2 is the list of possible pairings of 2 equations
  %--- pick one of these pairings
  l1:=nil;  
  %--- l1 is the list of equations which still can be reduced
  %--- and where f=car get(equ,'allvarfcts), i.e. equations
  %--- which must not be used for generating new equations
  %
  %--- each l in l2 has the form
  %--- (((e4 (x 2 y)) (e5 (z)) (df f x 2 y 2 z) nil) . num_value)
  for each l in l2 do <<
    f:=caddar l;
    if pairp f then f:=cadr f;
    if (caaar l = cadddr car l) and  % if caaar  l will be subst.
       get(caaar l,'allvarfcts) and
       (f=car get(caaar l,'allvarfcts)) 
    then l1:=union(list(caaar l),l1);
    if (caadar l = cadddr car l) and % if caadar l will be subst.
       get(caadar l,'allvarfcts) and
       (f=car get(caadar l,'allvarfcts))
    then l1:=union(list(caadar l),l1);
  >>;
  %--- Test that no new equation will be generated from an
  %--- equation from which the leading derivative can still be
  %--- reduced
  for each l in l2 do
  if ((cadaar l = nil)                              or 
      (cadr cadar l = nil)                          or
      (freeof(l1,caaar  l) and freeof(l1,caadar l))    ) and
     dec_info_leq(l,min) 
  then min:=l;
  if min then return dec_put_info(car min,rl)$
end$

symbolic procedure elimin(p,q,x)$
begin
 scalar dgs,s,flg,quoti,lts,
        ltp,dgp,fcpp,fcqp,fcsp,
        ltq,dgq,fcpq,fcqq,fcsq$

 ltp:=reval{'LTERM,p,x};
 ltq:=reval{'LTERM,q,x};
 dgp:=reval{'DEG,ltp,x};
 dgq:=reval{'DEG,ltq,x};

 if dgp > dgq then << flg:=t; dgs:=dgq >>
              else dgs:=dgp$

 fcpp:=1; fcpq:=0;
 fcqq:=1; fcqp:=0;

 while dgs neq 0 do <<
  quoti:=reval{'QUOTIENT,ltp,ltq};
  s:=reval{'PLUS,{'TIMES,p,{'DEN,quoti}},
                 {'MINUS,{'TIMES,q,{'NUM,quoti}}}}$
  lts:=reval{'LTERM,s,x}$
  dgs:=reval{'DEG,lts,x}$
  fcsp:=reval{'PLUS,{'TIMES,fcpp,{'DEN,quoti}},
                    {'MINUS,{'TIMES,fcqp,{'NUM,quoti}}}}$
  fcsq:=reval{'PLUS,{'TIMES,fcpq,{'DEN,quoti}},
                    {'MINUS,{'TIMES,fcqq,{'NUM,quoti}}}}$

  if flg=t then <<
   p:=s;
   ltp:=lts;
   dgp:=dgs;
   fcpp:=fcsp;
   fcpq:=fcsq;
   if dgq>dgp then flg:=nil
  >>       else <<
   q:=s;
   ltq:=lts;
   dgq:=dgs;
   fcqp:=fcsp;
   fcqq:=fcsq;
   if dgp>dgq then flg:=t
  >>$
 >>;
 quoti:=reval{'GCD,fcsp,fcsq};
 return {reval{'QUOTIENT,s   ,quoti},
         reval{'QUOTIENT,fcsp,quoti},
         reval{'QUOTIENT,fcsq,quoti} }
end$

symbolic procedure dec_reduction(l,ftem,vl,rl)$
% do a reduction step or a cross differentiation either
% l has the form  ((df f x 2 y 2 z) e4 e5 nil)
% if rl then one pde must be simplified with the help of 
% another one and reduce its length
begin scalar p,q,ld,a,s,ip,iq,f,dwsa,dwla,dwlb,el,h$
  ld:=car l$    % ld=deriv. to be elim. after differentiation
  l:=cdr l$     % l=list of equations
  f:=if pairp ld then cadr ld
                 else ld$
  p:=get(car l,'val)$
  q:=get(cadr l,'val)$
  ip:=get(car l,'dec_info)$
  iq:=get(cadr l,'dec_info)$
  s:=caddr l;               % s is the equation to be deleted
  l:=list(car l,cadr l);
  if print_ and ((null rl and tr_decouple) or 
                 (     rl and tr_redlength)  ) then
  <<terpri()$write "  first pde ",car l,": "$
    typeeq car l$ 
    if ip then write "is diff. wrt. ",ip,","
          else write "is not differentiated,"$
    write "  second pde ",cadr l,": "$
    typeeq cadr l$
    if iq then write "is diff. wrt. ",iq," "
          else write "is not differentiated, "$
    write"to eliminate ",mathprint ld$
  >>$
  if iq then q:=simplifypde(reval cons('DF,cons(q,iq)),ftem,nil)$
  if ip then p:=simplifypde(reval cons('DF,cons(p,ip)),ftem,nil)$
%  h:=reval !*q2a simpresultant list(p,q,ld)$
  h:=car elimin(p,q,ld)$
  if h and not zerop h 
  then a:=mkeq(h,ftem,vl,allflags_,t,list(0))$
  % Orderings ! Problem: we need to understand the internals of this
  % functions to a better degree.

  if a and rl then <<
    if get(s,'length)*length_inc<get(a,'length) then <<
      add_both_dec_with(f,car l,cadr l,t);
      if tr_redlength and print_ then 
      <<%terpri()$
        write "  The resulting pde "$typeeq a$
        write"  would have had ",get(a,'length),
        " terms compared to ",get(s,'length)," terms">>;
      setprop(a,nil)$
      nequ_:=sub1 nequ_$
      return nil
    >>
  >>$

  if a and member(get(a,'val),dec_hist_list) then <<
    setprop(a,nil)$
    nequ_:=sub1 nequ_$
    add_both_dec_with(f,car l,cadr l,rl);
    if print_ and ((rl and tr_redlength) or (null rl and tr_decouple)) 
    then <<terpri()$write "the resulting pde would lead to a cycle"$>> >>
  else <<
    if print_ then <<
%      terpri()$write if rl then "Length-reduction" else "Decoupling",
%                     " step with "$
%      fctprint l$write" : "$
%      terpri()$write if rl then "Length-reduction:" else "Decoupling:"$
      write" Eliminate "$
      mathprint ld$
      write"from ",if ip then cons('DF,cons(car l,ip))
                         else car  l, " and ",
                   if iq then cons('DF,cons(cadr l,iq))
                         else cadr l, ". "$
      if a then <<
        if s then << 
          write s,": "$
          typeeq s; %terpri()$
          write"is replaced by ",a,": "
        >>   else write a," is added: "$
        typeeq a
      >>     else if s then write s," is deleted"
    >>$
    if s then <<      % reduction 
      l:=delete(s,l)$
      if not (ip or iq) then << 
        % The equations wrt which s has already been decoupled
        % are to be listed under dec_with wrt to the equation
        % of both that is kept which is car l
        % purpose: These decouplings should not be done again
        % with car l as this may result in a cycle
        dwsa:=get(    s,'dec_with)$  
        dwla:=get(car l,'dec_with)$  
        for each el in dwsa do <<           
          dwlb:=assoc(car el,dwla)$
          dwla:=delete(dwlb,dwla)$
          if dwlb then dwlb:=cons(car el,union(cdr el,cdr dwlb))
                  else dwlb:=el$
          dwla:=cons(dwlb,dwla)
        >>$ 
        put(car l,'dec_with,dwla)$
      >>$
      % If other than the leading derivatives are reduced then
      % the new equ. a must inherit 'dec_with from equ. s
      if a and get(a,'derivs) and 
         (car get(a,'derivs) = car get(s,'derivs)) then <<
        dwsa:=get(s,'dec_with)$
        put(a,'dec_with,dwsa)
      >>$
      if dec_hist>0 then <<
        if length dec_hist_list>dec_hist then
           dec_hist_list:=cdr dec_hist_list$
        dec_hist_list:=reverse cons(get(s,'val),reverse dec_hist_list)$
      >>$
      setprop(s,nil)>>$
    if a then l:=append(l,list a)$
  >>$
  return l$
end$

symbolic procedure dec_fct_check(a,l)$
% checks, if a function occurs in only one pde
begin scalar ft$
  ft:=get(a,'fcts)$
  while ft and l do
    <<if flagp(car l,'to_decoup) then
        ft:=setdiff(ft,get(car l,'fcts))$
    l:=cdr l>>$
  if ft then remflag1(a,'to_decoup)$
return ft$
end$

symbolic procedure dec_one_step(l,ftem,vl,rl)$
% do one decouple step for f; i.e. select 2 pdes, differentiate them
% and add the new pde or replace an original one
% if rl then one pde must be simplified with the help of 
% another one and reduce its length
begin scalar l1,l2$
  if not expert_mode then l1:=l
                     else <<
    l1:=selectpdes(l,2)$
    drop_dec_with(car  l1,cadr l1,rl)$
    drop_dec_with(cadr l1,car  l1,rl)$
  >>$
again:
  l1:=dec_and_fct_select(l1,vl,rl)$
  if not l1 then l:=nil
  else 
    if (l2:=dec_reduction(l1,ftem,vl,rl)) then
      <<for each a in cdr l1 do l:=delete(a,l)$
        for each a in l2 do l:=eqinsert(a,l)$
        for each a in l2 do dec_fct_check(a,l)$
      >>
    else if rl then <<l1:=l;goto again>>
               else l:=nil$
  return l$
end$

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%  procedures for decoupling of similar pde                        %
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

%symbolic procedure rel_length_diff(p,q)$
%% print length difference in pro cent
%(abs(get(p,'length)-get(q,'length))*100)/
%   (get(p,'length)+get(q,'length))$

%symbolic procedure nearly_same(p,q)$
%begin scalar lp,lq$
% lp:=get(p,'fcts)$
% lq:=get(q,'fcts)$ 
% if null setdiff(get(p,'allvarfcts),get(q,'allvarfcts)) and
%    null setdiff(get(q,'allvarfcts),get(p,'allvarfcts)) and
%    ((length setdiff(lp,lq)+length setdiff(lq,lp))*100<
%    (length lp+length lq)*same_fcts) then
%    <<lp:=get(p,'derivs)$
%      lq:=get(q,'derivs)$
%      if (length setdiff(lp,lq)+length setdiff(lq,lp))*100<
%         (length lp+length lq)*same_derivs then return t>>$
%end$
       
%symbolic procedure get_same_pdes(pdes)$
%begin scalar l,n,res$
% while pdes do
%  <<l:=cdr pdes$
%  while l do
%    if (n:=rel_length_diff(car pdes,car l))<=same_length then
%       if nearly_same(car pdes,car l) then
%          <<res:=list(car pdes,car l)$
%            l:=nil>>
%       else l:=cdr l
%    else if n>5*same_length then l:=nil
%                            else l:=cdr l$
%  if res then pdes:=nil
%         else pdes:=cdr pdes
%  >>$
% return res$
%end$ 

endmodule$
end$
