element.red 2.5 KB

123456789101112131415161718192021222324252627282930313233343536373839404142434445464748495051525354555657585960616263646566676869707172737475767778798081
  1. module element;
  2. % Generate a random integral element
  3. % Author: David Hartley
  4. Comment. At present, the Cartan-Kaehler construction is used, as by
  5. Wahlquist, to reduce the problem to linear algebra. This fails for
  6. non-involutive systems.
  7. endcomment;
  8. put('integral_element,'rtypefn,'quotelist);
  9. put('integral_element,'listfn,'intelteval);
  10. symbolic procedure intelteval(u,v);
  11. % u:{eds}, v:bool -> intelteval:list of prefix
  12. if length u neq 1 then
  13. rerror(eds,000,"Wrong number of arguments to integral_element")
  14. else if not edsp(u := reval car u) then typerr(u,"EDS")
  15. else !*sys2a1(edscall intelt u,v);
  16. symbolic procedure intelt s;
  17. % s:eds -> intelt:sys
  18. % Produce an arbitrary integral element of s using the Cartan-Kaehler
  19. % construction.
  20. begin scalar g,v,a,h,z;
  21. s := closure s;
  22. g := gbsys s;
  23. % reduction in next lines ok since lpows g = prlkrns s
  24. v := foreach f in nonpfaffpart eds_sys s join
  25. if f := xreduce(f,eds_sys g) then {f};
  26. % get polar systems
  27. h := reversip foreach w on reverse indkrns s collect
  28. foreach f in v join
  29. foreach c in ordcomb(cdr w,degreepf f - 1) join
  30. if c := xcoeff(f,car w . c) then {lc c};
  31. % get graded variable list
  32. a := v := {};
  33. foreach w in indkrns s do
  34. << v := setdiff(foreach f in eds_sys g collect
  35. mvar numr lc xcoeff(f,{w}),a) . v;
  36. a := append(car v,a) >>;
  37. v := reverse v;
  38. % solve polar systems
  39. foreach x in pair(h,v) do
  40. << v := cdr x;
  41. x := foreach f in car x join if numr(f := subsq(f,z)) then {f};
  42. edsdebug("Polar system",x,'sq);
  43. z := append(edsransolve(x,v),z) >>;
  44. return foreach f in eds_sys g collect pullbackpf(f,z);
  45. end;
  46. symbolic procedure edsransolve(x,v);
  47. % x:list of sq, v:list of kernel -> edsransolve:map
  48. begin
  49. x := edssolve(x,v);
  50. if null x then
  51. rerror(eds,000,"Singular system in integral_element");
  52. if length x > 1 or null caar x then
  53. rerror(eds,000,"Bad system in integral_element");
  54. x := car cdr car x; % get the map part of first solution
  55. v := setdiff(v,foreach m in x collect car m);
  56. edsverbose({length v,"free variables"},nil,nil);
  57. v := foreach c in v collect c . sparserandom 5;
  58. x := nconc(pullbackmap(x,v),v);
  59. edsdebug("Solution",x,'map);
  60. return x;
  61. end;
  62. symbolic procedure sparserandom n;
  63. if random 100 < 0 then 0 else random(2*n+1)-n;
  64. endmodule;
  65. end;