wneuper/isac: src/sml/IsacKnowledge/Rational.ML@c58b3d2ee7a9

     1 (*.calculate in rationals: gcd, lcm, etc.

     2    (c) Stefan Karnel 2002

     3    Institute for Mathematics D and Institute for Software Technology,

     4    TU-Graz SS 2002

     5    Use is subject to license terms.

     7 use"IsacKnowledge/Rational.ML";

     8 use"Rational.ML";

    10 remove_thy"Rational";

    11 use_thy"IsacKnowledge/Isac";

    12 ****************************************************************.*)

    14 (*.*****************************************************************

    15   Remark on notions in the documentation below:

    16     referring to the remark on 'polynomials' in Poly.sml we use

    17     [2] 'polynomial' normalform (Polynom)

    18     [3] 'expanded_term' normalform (Ausmultiplizierter Term),

    19     where normalform [2] is a special case of [3], i.e. [3] implies [2].

    20     Instead of

    21       'fraction with numerator and nominator both in normalform [2]'

    22       'fraction with numerator and nominator both in normalform [3]'

    23     we say:

    24       'fraction in normalform [2]'

    25       'fraction in normalform [3]'

    26     or

    27       'fraction [2]'

    28       'fraction [3]'.

    29     a 'simple fraction' is a term with '/' as outmost operator and

    30     numerator and nominator in normalform [2] or [3].

    31 ****************************************************************.*)

    33 signature RATIONALI =

    34 sig

    35   type mv_monom

    36   type mv_poly

    37   val add_fraction_ : theory -> term -> (term * term list) option

    38   val add_fraction_p_ : theory -> term -> (term * term list) option

    39   val calculate_Rational : rls

    40   val calc_rat_erls:rls

    41   val cancel : rls

    42   val cancel_ : theory -> term -> (term * term list) option

    43   val cancel_p : rls

    44   val cancel_p_ : theory -> term -> (term * term list) option

    45   val common_nominator : rls

    46   val common_nominator_ : theory -> term -> (term * term list) option

    47   val common_nominator_p : rls

    48   val common_nominator_p_ : theory -> term -> (term * term list) option

    49   val eval_is_expanded : string -> 'a -> term -> theory ->

    50 			 (string * term) option

    51   val expanded2polynomial : term -> term option

    52   val factout_ : theory -> term -> (term * term list) option

    53   val factout_p_ : theory -> term -> (term * term list) option

    54   val is_expanded : term -> bool

    55   val is_polynomial : term -> bool

    57   val mv_gcd : (int * int list) list -> mv_poly -> mv_poly

    58   val mv_lcm : mv_poly -> mv_poly -> mv_poly

    60   val norm_expanded_rat_ : theory -> term -> (term * term list) option

    61 (*WN0602.2.6.pull into struct !!!

    62   val norm_Rational : rls(*.normalizes an arbitrary rational term without

    63                            roots into a simple and canceled fraction

    64                            with normalform [2].*)

    65 *)

    66 (*val norm_rational_p : 19.10.02 missing FIXXXXXXXXXXXXME

    67       rls               (*.normalizes an rational term [2] without

    68                            roots into a simple and canceled fraction

    69                            with normalform [2].*)

    70 *)

    71   val norm_rational_ : theory -> term -> (term * term list) option

    72   val polynomial2expanded : term -> term option

    73   val rational_erls :

    74       rls             (*.evaluates an arbitrary rational term with numerals.*)

    76 (*WN0210???SK: fehlen Funktionen, die exportiert werden sollen ? *)

    77 end

    79 (*.**************************************************************************

    80 survey on the functions

    81 ~~~~~~~~~~~~~~~~~~~~~~~

    82  [2] 'polynomial'   :rls               | [3]'expanded_term':rls

    83 --------------------:------------------+-------------------:-----------------

    84  factout_p_         :                  | factout_          :

    85  cancel_p_          :                  | cancel_           :

    86                     :cancel_p          |                   :cancel

    87 --------------------:------------------+-------------------:-----------------

    88  common_nominator_p_:                  | common_nominator_ :

    89                     :common_nominator_p|                   :common_nominator

    90  add_fraction_p_    :                  | add_fraction_     :

    91 --------------------:------------------+-------------------:-----------------

    92 ???SK                 :norm_rational_p   |                   :norm_rational

    94 This survey shows only the principal functions for reuse, and the identifiers

    95 of the rls exported. The list below shows some more useful functions.

    98 conversion from Isabelle-term to internal representation

    99 ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~

   101 ... BITTE FORTSETZEN ...

   103 polynomial2expanded = ...

   104 expanded2polynomial = ...

   106 remark: polynomial2expanded o expanded2polynomial = I,

   107         where 'o' is function chaining, and 'I' is identity WN0210???SK

   109 functions for greatest common divisor and canceling

   110 ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~

   111 mv_gcd

   112 factout_

   113 factout_p_

   114 cancel_

   115 cancel_p_

   117 functions for least common multiple and addition of fractions

   118 ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~

   119 mv_lcm

   120 common_nominator_

   121 common_nominator_p_

   122 add_fraction_       (*.add 2 or more fractions.*)

   123 add_fraction_p_     (*.add 2 or more fractions.*)

   125 functions for normalform of rationals

   126 ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~

   127 WN0210???SK interne Funktionen f"ur norm_rational:

   128           schaffen diese SML-Funktionen wirklich ganz allgemeine Terme

   130 norm_rational_

   131 norm_expanded_rat_

   133 **************************************************************************.*)

   136 (*##*)

   137 structure RationalI : RATIONALI =

   138 struct

   139 (*##*)

   141 (*. gcd of integers .*)

   142 (* die gcd Funktion von Isabelle funktioniert nicht richtig !!! *)

   143 fun gcd_int a b = if b=0 then a

   144 		  else  gcd_int b (a mod b);

   146 (*. univariate polynomials (uv) .*)

   147 (*. univariate polynomials are represented as a list of the coefficent in reverse maximum degree order .*)

   148 (*. 5 * x^5 + 4 * x^3 + 2 * x^2 + x + 19 => [19,1,2,4,0,5] .*)

   149 type uv_poly = int list;

   151 (*. adds two uv polynomials .*)

   152 fun uv_mod_add_poly ([]:uv_poly,p2:uv_poly) = p2:uv_poly

   153   | uv_mod_add_poly (p1,[]) = p1

   154   | uv_mod_add_poly (x::p1,y::p2) = (x+y)::(uv_mod_add_poly(p1,p2));

   156 (*. multiplies a uv polynomial with a skalar s .*)

   157 fun uv_mod_smul_poly ([]:uv_poly,s:int) = []:uv_poly

   158   | uv_mod_smul_poly (x::p,s) = (x*s)::(uv_mod_smul_poly(p,s));

   160 (*. calculates the remainder of a polynomial divided by a skalar s .*)

   161 fun uv_mod_rem_poly ([]:uv_poly,s) = []:uv_poly

   162   | uv_mod_rem_poly (x::p,s) = (x mod s)::(uv_mod_smul_poly(p,s));

   164 (*. calculates the degree of a uv polynomial .*)

   165 fun uv_mod_deg ([]:uv_poly) = 0

   166   | uv_mod_deg p = length(p)-1;

   168 (*. calculates the remainder of x/p and represents it as value between -p/2 and p/2 .*)

   169 fun uv_mod_mod2(x,p)=

   170     let

   171 	val y=(x mod p);

   172     in

   173 	if (y)>(p div 2) then (y)-p else

   174 	    (

   175 	     if (y)<(~p div 2) then p+(y) else (y)

   176 	     )

   177     end;

   179 (*.calculates the remainder for each element of a integer list divided by p.*)

   180 fun uv_mod_list_modp [] p = []

   181   | uv_mod_list_modp (x::xs) p = (uv_mod_mod2(x,p))::(uv_mod_list_modp xs p);

   183 (*. appends an integer at the end of a integer list .*)

   184 fun uv_mod_null (p1:int list,0) = p1

   185   | uv_mod_null (p1:int list,n1:int) = uv_mod_null(p1,n1-1) @ [0];

   187 (*. uv polynomial division, result is (quotient, remainder) .*)

   188 (*. only for uv_mod_divides .*)

   189 (* FIXME: Division von x^9+x^5+1 durch x-1000 funktioniert nicht integer zu klein  *)

   190 fun uv_mod_pdiv (p1:uv_poly) ([]:uv_poly) = raise error ("RATIONALS_UV_MOD_PDIV_EXCEPTION: division by zero")

   191   | uv_mod_pdiv p1 [x] =

   192     let

   193 	val xs=ref [];

   194     in

   195 	if x<>0 then

   196 	    (

   197 	     xs:=(uv_mod_rem_poly(p1,x));

   198 	     while length(!xs)>0 andalso hd(!xs)=0 do xs:=tl(!xs)

   199 	     )

   200 	else raise error ("RATIONALS_UV_MOD_PDIV_EXCEPTION: division by zero");

   201 	([]:uv_poly,!xs:uv_poly)

   202     end

   203   | uv_mod_pdiv p1 p2 =

   204     let

   205 	val n= uv_mod_deg(p2);

   206 	val m= ref (uv_mod_deg(p1));

   207 	val p1'=ref (rev(p1));

   208 	val p2'=(rev(p2));

   209 	val lc2=hd(p2');

   210 	val q=ref [];

   211 	val c=ref 0;

   212 	val output=ref ([],[]);

   213     in

   214 	(

   215 	 if (!m)=0 orelse p2=[0] then raise error ("RATIONALS_UV_MOD_PDIV_EXCEPTION: Division by zero")

   216 	 else

   217 	     (

   218 	      if (!m)<n then

   219 		  (

   220 		   output:=([0],p1)

   221 		   )

   222 	      else

   223 		  (

   224 		   while (!m)>=n do

   225 		       (

   226 			c:=hd(!p1') div hd(p2');

   227 			if !c<>0 then

   228 			    (

   229 			     p1':=uv_mod_add_poly(!p1',uv_mod_null(uv_mod_smul_poly(p2',~(!c)),!m-n));

   230 			     while length(!p1')>0 andalso hd(!p1')=0  do p1':= tl(!p1');

   231 			     m:=uv_mod_deg(!p1')

   232 			     )

   233 			else m:=0

   234 			);

   235     		   output:=(rev(!q),rev(!p1'))

   236 		   )

   237 	      );

   238 	     !output

   239 	 )

   240     end;

   242 (*. divides p1 by p2 in Zp .*)

   243 fun uv_mod_pdivp (p1:uv_poly) (p2:uv_poly) p =

   244     let

   245 	val n=uv_mod_deg(p2);

   246 	val m=ref (uv_mod_deg(uv_mod_list_modp p1 p));

   247 	val p1'=ref (rev(p1));

   248 	val p2'=(rev(uv_mod_list_modp p2 p));

   249 	val lc2=hd(p2');

   250 	val q=ref [];

   251 	val c=ref 0;

   252 	val output=ref ([],[]);

   253     in

   254 	(

   255 	 if (!m)=0 orelse p2=[0] then raise error ("RATIONALS_UV_MOD_PDIVP_EXCEPTION: Division by zero")

   256 	 else

   257 	     (

   258 	      if (!m)<n then

   259 		  (

   260 		   output:=([0],p1)

   261 		   )

   262 	      else

   263 		  (

   264 		   while (!m)>=n do

   265 		       (

   266 			c:=uv_mod_mod2(hd(!p1')*(power lc2 1), p);

   267 			q:=(!c)::(!q);

   268 			p1':=uv_mod_list_modp(tl(uv_mod_add_poly(uv_mod_smul_poly(!p1',lc2),

   269 								  uv_mod_smul_poly(uv_mod_smul_poly(p2',hd(!p1')),~1)))) p;

   270 			m:=(!m)-1

   271 			);

   273 		   while !p1'<>[] andalso hd(!p1')=0 do

   274 		       (

   275 			p1':=tl(!p1')

   276 			);

   278     		   output:=(rev(uv_mod_list_modp (!q) (p)),rev(!p1'))

   279 		   )

   280 	      );

   281 	     !output:uv_poly * uv_poly

   282 	 )

   283     end;

   285 (*. calculates the remainder of p1/p2 .*)

   286 fun uv_mod_prest (p1:uv_poly) ([]:uv_poly) = raise error("UV_MOD_PREST_EXCEPTION: Division by zero")

   287   | uv_mod_prest [] p2 = []:uv_poly

   288   | uv_mod_prest p1 p2 = (#2(uv_mod_pdiv p1 p2));

   290 (*. calculates the remainder of p1/p2 in Zp .*)

   291 fun uv_mod_prestp (p1:uv_poly) ([]:uv_poly) p= raise error("UV_MOD_PRESTP_EXCEPTION: Division by zero")

   292   | uv_mod_prestp [] p2 p= []:uv_poly

   293   | uv_mod_prestp p1 p2 p = #2(uv_mod_pdivp p1 p2 p);

   295 (*. calculates the content of a uv polynomial .*)

   296 fun uv_mod_cont ([]:uv_poly) = 0

   297   | uv_mod_cont (x::p)= gcd_int x (uv_mod_cont(p));

   299 (*. divides each coefficient of a uv polynomial by y .*)

   300 fun uv_mod_div_list (p:uv_poly,0) = raise error("UV_MOD_DIV_LIST_EXCEPTION: Division by zero")

   301   | uv_mod_div_list ([],y)   = []:uv_poly

   302   | uv_mod_div_list (x::p,y) = (x div y)::uv_mod_div_list(p,y);

   304 (*. calculates the primitiv part of a uv polynomial .*)

   305 fun uv_mod_pp ([]:uv_poly) = []:uv_poly

   306   | uv_mod_pp p =

   307     let

   308 	val c=ref 0;

   309     in

   310 	(

   311 	 c:=uv_mod_cont(p);

   313 	 if !c=0 then raise error ("RATIONALS_UV_MOD_PP_EXCEPTION: content is 0")

   314 	 else uv_mod_div_list(p,!c)

   315 	)

   316     end;

   318 (*. gets the leading coefficient of a uv polynomial .*)

   319 fun uv_mod_lc ([]:uv_poly) = 0

   320   | uv_mod_lc p  = hd(rev(p));

   322 (*. calculates the euklidean polynomial remainder sequence in Zp .*)

   323 fun uv_mod_prs_euklid_p(p1:uv_poly,p2:uv_poly,p)=

   324     let

   325 	val f =ref [];

   326 	val f'=ref p2;

   327 	val fi=ref [];

   328     in

   329 	(

   330 	 f:=p2::p1::[];

   331  	 while uv_mod_deg(!f')>0 do

   332 	     (

   333 	      f':=uv_mod_prestp (hd(tl(!f))) (hd(!f)) p;

   334 	      if (!f')<>[] then

   335 		  (

   336 		   fi:=(!f');

   337 		   f:=(!fi)::(!f)

   338 		   )

   339 	      else ()

   340 	      );

   341 	      (!f)

   343 	 )

   344     end;

   346 (*. calculates the gcd of p1 and p2 in Zp .*)

   347 fun uv_mod_gcd_modp ([]:uv_poly) (p2:uv_poly) p = p2:uv_poly

   348   | uv_mod_gcd_modp p1 [] p= p1

   349   | uv_mod_gcd_modp p1 p2 p=

   350     let

   351 	val p1'=ref[];

   352 	val p2'=ref[];

   353 	val pc=ref[];

   354 	val g=ref [];

   355 	val d=ref 0;

   356 	val prs=ref [];

   357     in

   358 	(

   359 	 if uv_mod_deg(p1)>=uv_mod_deg(p2) then

   360 	     (

   361 	      p1':=uv_mod_list_modp (uv_mod_pp(p1)) p;

   362 	      p2':=uv_mod_list_modp (uv_mod_pp(p2)) p

   363 	      )

   364 	 else

   365 	     (

   366 	      p1':=uv_mod_list_modp (uv_mod_pp(p2)) p;

   367 	      p2':=uv_mod_list_modp (uv_mod_pp(p1)) p

   368 	      );

   369 	 d:=uv_mod_mod2((gcd_int (uv_mod_cont(p1))) (uv_mod_cont(p2)), p) ;

   370 	 if !d>(p div 2) then d:=(!d)-p else ();

   372 	 prs:=uv_mod_prs_euklid_p(!p1',!p2',p);

   374 	 if hd(!prs)=[] then pc:=hd(tl(!prs))

   375 	 else pc:=hd(!prs);

   377 	 g:=uv_mod_smul_poly(uv_mod_pp(!pc),!d);

   378 	 !g

   379 	 )

   380     end;

   382 (*. calculates the minimum of two real values x and y .*)

   383 fun uv_mod_r_min(x,y):BasisLibrary.Real.real = if x>y then y else x;

   385 (*. calculates the minimum of two integer values x and y .*)

   386 fun uv_mod_min(x,y) = if x>y then y else x;

   388 (*. adds the squared values of a integer list .*)

   389 fun uv_mod_add_qu [] = 0.0

   390   | uv_mod_add_qu (x::p) =  BasisLibrary.Real.fromInt(x)*BasisLibrary.Real.fromInt(x) + uv_mod_add_qu p;

   392 (*. calculates the euklidean norm .*)

   393 fun uv_mod_norm ([]:uv_poly) = 0.0

   394   | uv_mod_norm p = Math.sqrt(uv_mod_add_qu(p));

   396 (*. multipies two values a and b .*)

   397 fun uv_mod_multi a b = a * b;

   399 (*. decides if x is a prim, the list contains all primes which are lower then x .*)

   400 fun uv_mod_prim(x,[])= false

   401   | uv_mod_prim(x,[y])=if ((x mod y) <> 0) then true

   402 		else false

   403   | uv_mod_prim(x,y::ys) = if uv_mod_prim(x,[y])

   404 			then

   405 			    if uv_mod_prim(x,ys) then true

   406 			    else false

   407 		    else false;

   409 (*. gets the first prime, which is greater than p and does not divide g .*)

   410 fun uv_mod_nextprime(g,p)=

   411     let

   412 	val list=ref [2];

   413 	val exit=ref 0;

   414 	val i = ref 2

   415     in

   416 	while (!i<p) do (* calculates the primes lower then p *)

   417 	    (

   418 	     if uv_mod_prim(!i,!list) then

   419 		 (

   420 		  if (p mod !i <> 0)

   421 		      then

   422 			  (

   423 			   list:= (!i)::(!list);

   424 			   i:= (!i)+1

   425 			   )

   426 		  else i:=(!i)+1

   427 		  )

   428 	     else i:= (!i)+1

   429 		 );

   430 	    i:=(p+1);

   431 	    while (!exit=0) do   (* calculate next prime which does not divide g *)

   432 	    (

   433 	     if uv_mod_prim(!i,!list) then

   434 		 (

   435 		  if (g mod !i <> 0)

   436 		      then

   437 			  (

   438 			   list:= (!i)::(!list);

   439 			   exit:= (!i)

   440 			   )

   441 		  else i:=(!i)+1

   442 		      )

   443 	     else i:= (!i)+1

   444 		 );

   445 	    !exit

   446     end;

   448 (*. decides if p1 is a factor of p2 in Zp .*)

   449 fun uv_mod_dividesp ([]:uv_poly) (p2:uv_poly) p= raise error("UV_MOD_DIVIDESP: Division by zero")

   450   | uv_mod_dividesp p1 p2 p= if uv_mod_prestp p2 p1 p = [] then true else false;

   452 (*. decides if p1 is a factor of p2 .*)

   453 fun uv_mod_divides ([]:uv_poly) (p2:uv_poly) = raise error("UV_MOD_DIVIDES: Division by zero")

   454   | uv_mod_divides p1 p2 = if uv_mod_prest p2 p1  = [] then true else false;

   456 (*. chinese remainder algorithm .*)

   457 fun uv_mod_cra2(r1,r2,m1,m2)=

   458     let

   459 	val c=ref 0;

   460 	val r1'=ref 0;

   461 	val d=ref 0;

   462 	val a=ref 0;

   463     in

   464 	(

   465 	 while (uv_mod_mod2((!c)*m1,m2))<>1 do

   466 	     (

   467 	      c:=(!c)+1

   468 	      );

   469 	 r1':= uv_mod_mod2(r1,m1);

   470 	 d:=uv_mod_mod2(((r2-(!r1'))*(!c)),m2);

   471 	 !r1'+(!d)*m1

   472 	 )

   473     end;

   475 (*. applies the chinese remainder algorithmen to the coefficients of x1 and x2 .*)

   476 fun uv_mod_cra_2 ([],[],m1,m2) = []

   477   | uv_mod_cra_2 ([],x2,m1,m2) = raise error("UV_MOD_CRA_2_EXCEPTION: invalid call x1")

   478   | uv_mod_cra_2 (x1,[],m1,m2) = raise error("UV_MOD_CRA_2_EXCEPTION: invalid call x2")

   479   | uv_mod_cra_2 (x1::x1s,x2::x2s,m1,m2) = (uv_mod_cra2(x1,x2,m1,m2))::(uv_mod_cra_2(x1s,x2s,m1,m2));

   481 (*. calculates the gcd of two uv polynomials p1' and p2' with the modular algorithm .*)

   482 fun uv_mod_gcd (p1':uv_poly) (p2':uv_poly) =

   483     let

   484 	val p1=ref (uv_mod_pp(p1'));

   485 	val p2=ref (uv_mod_pp(p2'));

   486 	val c=gcd_int (uv_mod_cont(p1')) (uv_mod_cont(p2'));

   487 	val temp=ref [];

   488 	val cp=ref [];

   489 	val qp=ref [];

   490 	val q=ref[];

   491 	val pn=ref 0;

   492 	val d=ref 0;

   493 	val g1=ref 0;

   494 	val p=ref 0;

   495 	val m=ref 0;

   496 	val exit=ref 0;

   497 	val i=ref 1;

   498     in

   499 	if length(!p1)>length(!p2) then ()

   500 	else

   501 	    (

   502 	     temp:= !p1;

   503 	     p1:= !p2;

   504 	     p2:= !temp

   505 	     );

   508 	d:=gcd_int (uv_mod_lc(!p1)) (uv_mod_lc(!p2));

   509 	g1:=uv_mod_lc(!p1)*uv_mod_lc(!p2);

   510 	p:=4;

   512 	m:=BasisLibrary.Real.ceil(2.0 *

   513 				  BasisLibrary.Real.fromInt(!d) *

   514 				  BasisLibrary.Real.fromInt(power 2 (uv_mod_min(uv_mod_deg(!p2),uv_mod_deg(!p1)))) *

   515 				  BasisLibrary.Real.fromInt(!d) *

   516 				  uv_mod_r_min(uv_mod_norm(!p1) / BasisLibrary.Real.fromInt(abs(uv_mod_lc(!p1))),

   517 					uv_mod_norm(!p2) / BasisLibrary.Real.fromInt(abs(uv_mod_lc(!p2)))));

   519 	while (!exit=0) do

   520 	    (

   521 	     p:=uv_mod_nextprime(!d,!p);

   522 	     cp:=(uv_mod_gcd_modp (uv_mod_list_modp(!p1) (!p)) (uv_mod_list_modp(!p2) (!p)) (!p)) ;

   523 	     if abs(uv_mod_lc(!cp))<>1 then  (* leading coefficient = 1 ? *)

   524 		 (

   525 		  i:=1;

   526 		  while (!i)<(!p) andalso (abs(uv_mod_mod2((uv_mod_lc(!cp)*(!i)),(!p)))<>1) do

   527 		      (

   528 		       i:=(!i)+1

   529 		       );

   530 		      cp:=uv_mod_list_modp (map (uv_mod_multi (!i)) (!cp)) (!p)

   531 		  )

   532 	     else ();

   534 	     qp:= ((map (uv_mod_multi (uv_mod_mod2(!d,!p)))) (!cp));

   536 	     if uv_mod_deg(!qp)=0 then (q:=[1]; exit:=1) else ();

   538 	     pn:=(!p);

   539 	     q:=(!qp);

   541 	     while !pn<= !m andalso !m>(!p) andalso !exit=0 do

   542 		 (

   543 		  p:=uv_mod_nextprime(!d,!p);

   544  		  cp:=(uv_mod_gcd_modp (uv_mod_list_modp(!p1) (!p)) (uv_mod_list_modp(!p2) (!p)) (!p));

   545  		  if uv_mod_lc(!cp)<>1 then  (* leading coefficient = 1 ? *)

   546  		      (

   547  		       i:=1;

   548  		       while (!i)<(!p) andalso ((uv_mod_mod2((uv_mod_lc(!q)*(!i)),(!p)))<>1) do

   549  			   (

   550  			    i:=(!i)+1

   551 		           );

   552 		       cp:=uv_mod_list_modp (map (uv_mod_multi (!i)) (!cp)) (!p)

   553  		      )

   554  		  else ();

   556 		  qp:=uv_mod_list_modp ((map (uv_mod_multi (uv_mod_mod2(!d,!p)))) (!cp)  ) (!p);

   557  		  if uv_mod_deg(!qp)>uv_mod_deg(!q) then

   558  		      (

   559  		       (*print("degree to high!!!\n")*)

   560  		       )

   561  		  else

   562  		      (

   563  		       if uv_mod_deg(!qp)=uv_mod_deg(!q) then

   564  			   (

   565  			    q:=uv_mod_cra_2(!q,!qp,!pn,!p);

   566 			    pn:=(!pn) * !p;

   567 			    q:=uv_mod_pp(uv_mod_list_modp (!q) (!pn)); (* found already gcd ? *)

   568 			    if (uv_mod_divides (!q) (p1')) andalso (uv_mod_divides (!q) (p2')) then (exit:=1) else ()

   569 		 	    )

   570 		       else

   571 			   (

   572 			    if  uv_mod_deg(!qp)<uv_mod_deg(!q) then

   573 				(

   574 				 pn:= !p;

   575 				 q:= !qp

   576 				 )

   577 			    else ()

   578 			    )

   579 		       )

   580 		  );

   581  	     q:=uv_mod_pp(uv_mod_list_modp (!q) (!pn));

   582 	     if (uv_mod_divides (!q) (p1')) andalso (uv_mod_divides (!q) (p2')) then exit:=1 else ()

   583 	     );

   584 	    uv_mod_smul_poly(!q,c):uv_poly

   585     end;

   587 (*. multivariate polynomials .*)

   588 (*. multivariate polynomials are represented as a list of the pairs,

   589  first is the coefficent and the second is a list of the exponents .*)

   590 (*. 5 * x^5 * y^3 + 4 * x^3 * z^2 + 2 * x^2 * y * z^3 - z - 19

   591  => [(5,[5,3,0]),(4,[3,0,2]),(2,[2,1,3]),(~1,[0,0,1]),(~19,[0,0,0])] .*)

   593 (*. global variables .*)

   594 (*. order indicators .*)

   595 val LEX_=0; (* lexicographical term order *)

   596 val GGO_=1; (* greatest degree order *)

   598 (*. datatypes for internal representation.*)

   599 type mv_monom = (int *      (*.coefficient or the monom.*)

   600 		 int list); (*.list of exponents)      .*)

   601 fun mv_monom2str (i, is) = "("^ int2str i^"," ^ ints2str' is ^ ")";

   603 type mv_poly = mv_monom list;

   604 fun mv_poly2str p = (strs2str' o (map mv_monom2str)) p;

   606 (*. help function for monom_greater and geq .*)

   607 fun mv_mg_hlp([]) = EQUAL

   608   | mv_mg_hlp(x::list)=if x<0 then LESS

   609 		    else if x>0 then GREATER

   610 			 else mv_mg_hlp(list);

   612 (*. adds a list of values .*)

   613 fun mv_addlist([]) = 0

   614   | mv_addlist(p1) = hd(p1)+mv_addlist(tl(p1));

   616 (*. tests if the monomial M1 is greater as the monomial M2 and returns a boolean value .*)

   617 (*. 2 orders are implemented LEX_/GGO_ (lexigraphical/greatest degree order) .*)

   618 fun mv_monom_greater((M1x,M1l):mv_monom,(M2x,M2l):mv_monom,order)=

   619     if order=LEX_ then

   620 	(

   621 	 if length(M1l)<>length(M2l) then raise error ("RATIONALS_MV_MONOM_GREATER_EXCEPTION: Order error")

   622 	 else if (mv_mg_hlp((map op- (M1l~~M2l)))<>GREATER) then false else true

   623 	     )

   624     else

   625 	if order=GGO_ then

   626 	    (

   627 	     if length(M1l)<>length(M2l) then raise error ("RATIONALS_MV_MONOM_GREATER_EXCEPTION: Order error")

   628 	     else

   629 		 if mv_addlist(M1l)=mv_addlist(M2l)  then if (mv_mg_hlp((map op- (M1l~~M2l)))<>GREATER) then false else true

   630 		 else if mv_addlist(M1l)>mv_addlist(M2l) then true else false

   631 	     )

   632 	else raise error ("RATIONALS_MV_MONOM_GREATER_EXCEPTION: Wrong Order");

   634 (*. tests if the monomial X is greater as the monomial Y and returns a order value (GREATER,EQUAL,LESS) .*)

   635 (*. 2 orders are implemented LEX_/GGO_ (lexigraphical/greatest degree order) .*)

   636 fun mv_geq order ((x1,x):mv_monom,(x2,y):mv_monom) =

   637 let

   638     val temp=ref EQUAL;

   639 in

   640     if order=LEX_ then

   641 	(

   642 	 if length(x)<>length(y) then

   643 	     raise error ("RATIONALS_MV_GEQ_EXCEPTION: Order error")

   644 	 else

   645 	     (

   646 	      temp:=mv_mg_hlp((map op- (x~~y)));

   647 	      if !temp=EQUAL then

   648 		  ( if x1=x2 then EQUAL

   649 		    else if x1>x2 then GREATER

   650 			 else LESS

   651 			     )

   652 	      else (!temp)

   653 	      )

   654 	     )

   655     else

   656 	if order=GGO_ then

   657 	    (

   658 	     if length(x)<>length(y) then

   659 	      raise error ("RATIONALS_MV_GEQ_EXCEPTION: Order error")

   660 	     else

   661 		 if mv_addlist(x)=mv_addlist(y) then

   662 		     (mv_mg_hlp((map op- (x~~y))))

   663 		 else if mv_addlist(x)>mv_addlist(y) then GREATER else LESS

   664 		     )

   665 	else raise error ("RATIONALS_MV_GEQ_EXCEPTION: Wrong Order")

   666 end;

   668 (*. cuts the first variable from a polynomial .*)

   669 fun mv_cut([]:mv_poly)=[]:mv_poly

   670   | mv_cut((x,[])::list) = raise error ("RATIONALS_MV_CUT_EXCEPTION: Invalid list ")

   671   | mv_cut((x,y::ys)::list)=(x,ys)::mv_cut(list);

   673 (*. leading power product .*)

   674 fun mv_lpp([]:mv_poly,order)  = []

   675   | mv_lpp([(x,y)],order) = y

   676   | mv_lpp(p1,order)  = #2(hd(rev(sort (mv_geq order) p1)));

   678 (*. leading monomial .*)

   679 fun mv_lm([]:mv_poly,order)  = (0,[]):mv_monom

   680   | mv_lm([x],order) = x

   681   | mv_lm(p1,order)  = hd(rev(sort (mv_geq order) p1));

   683 (*. leading coefficient in term order .*)

   684 fun mv_lc2([]:mv_poly,order)  = 0

   685   | mv_lc2([(x,y)],order) = x

   686   | mv_lc2(p1,order)  = #1(hd(rev(sort (mv_geq order) p1)));

   689 (*. reverse the coefficients in mv polynomial .*)

   690 fun mv_rev_to([]:mv_poly) = []:mv_poly

   691   | mv_rev_to((c,e)::xs) = (c,rev(e))::mv_rev_to(xs);

   693 (*. leading coefficient in reverse term order .*)

   694 fun mv_lc([]:mv_poly,order)  = []:mv_poly

   695   | mv_lc([(x,y)],order) = mv_rev_to(mv_cut(mv_rev_to([(x,y)])))

   696   | mv_lc(p1,order)  =

   697     let

   698 	val p1o=ref (rev(sort (mv_geq order) (mv_rev_to(p1))));

   699 	val lp=hd(#2(hd(!p1o)));

   700 	val lc=ref [];

   701     in

   702 	(

   703 	 while (length(!p1o)>0 andalso hd(#2(hd(!p1o)))=lp) do

   704 	     (

   705 	      lc:=hd(mv_cut([hd(!p1o)]))::(!lc);

   706 	      p1o:=tl(!p1o)

   707 	      );

   708 	 if !lc=[] then raise error ("RATIONALS_MV_LC_EXCEPTION: lc is empty") else ();

   709 	 mv_rev_to(!lc)

   710 	 )

   711     end;

   713 (*. compares two powerproducts .*)

   714 fun mv_monom_equal((_,xlist):mv_monom,(_,ylist):mv_monom) = (foldr and_) (((map op=) (xlist~~ylist)),true);

   716 (*. help function for mv_add .*)

   717 fun mv_madd([]:mv_poly,[]:mv_poly,order) = []:mv_poly

   718   | mv_madd([(0,_)],p2,order) = p2

   719   | mv_madd(p1,[(0,_)],order) = p1

   720   | mv_madd([],p2,order) = p2

   721   | mv_madd(p1,[],order) = p1

   722   | mv_madd(p1,p2,order) =

   723     (

   724      if mv_monom_greater(hd(p1),hd(p2),order)

   725 	 then hd(p1)::mv_madd(tl(p1),p2,order)

   726      else if mv_monom_equal(hd(p1),hd(p2))

   727 	      then if mv_lc2(p1,order)+mv_lc2(p2,order)<>0

   728 		       then (mv_lc2(p1,order)+mv_lc2(p2,order),mv_lpp(p1,order))::mv_madd(tl(p1),tl(p2),order)

   729 		   else mv_madd(tl(p1),tl(p2),order)

   730 	  else hd(p2)::mv_madd(p1,tl(p2),order)

   731 	      )

   733 (*. adds two multivariate polynomials .*)

   734 fun mv_add([]:mv_poly,p2:mv_poly,order) = p2

   735   | mv_add(p1,[],order) = p1

   736   | mv_add(p1,p2,order) = mv_madd(rev(sort (mv_geq order) p1),rev(sort (mv_geq order) p2), order);

   738 (*. monom multiplication .*)

   739 fun mv_mmul((x1,y1):mv_monom,(x2,y2):mv_monom)=(x1*x2,(map op+) (y1~~y2)):mv_monom;

   741 (*. deletes all monomials with coefficient 0 .*)

   742 fun mv_shorten([]:mv_poly,order) = []:mv_poly

   743   | mv_shorten(x::xs,order)=mv_madd([x],mv_shorten(xs,order),order);

   745 (*. zeros a list .*)

   746 fun mv_null2([])=[]

   747   | mv_null2(x::l)=0::mv_null2(l);

   749 (*. multiplies two multivariate polynomials .*)

   750 fun mv_mul([]:mv_poly,[]:mv_poly,_) = []:mv_poly

   751   | mv_mul([],y::p2,_) = [(0,mv_null2(#2(y)))]

   752   | mv_mul(x::p1,[],_) = [(0,mv_null2(#2(x)))]

   753   | mv_mul(x::p1,y::p2,order) = mv_shorten(rev(sort (mv_geq order) (mv_mmul(x,y) :: (mv_mul(p1,y::p2,order) @

   754 									    mv_mul([x],p2,order)))),order);

   756 (*. gets the maximum value of a list .*)

   757 fun mv_getmax([])=0

   758   | mv_getmax(x::p1)= let

   759 		       val m=mv_getmax(p1);

   760 		   in

   761 		       if m>x then m

   762 		       else x

   763 		   end;

   764 (*. calculates the maximum degree of an multivariate polynomial .*)

   765 fun mv_grad([]:mv_poly) = 0

   766   | mv_grad(p1:mv_poly)= mv_getmax((map mv_addlist) ((map #2) p1));

   768 (*. converts the sign of a value .*)

   769 fun mv_minus(x)=(~1) * x;

   771 (*. converts the sign of all coefficients of a polynomial .*)

   772 fun mv_minus2([]:mv_poly)=[]:mv_poly

   773   | mv_minus2(p1)=(mv_minus(#1(hd(p1))),#2(hd(p1)))::(mv_minus2(tl(p1)));

   775 (*. searches for a negativ value in a list .*)

   776 fun mv_is_negativ([])=false

   777   | mv_is_negativ(x::xs)=if x<0 then true else mv_is_negativ(xs);

   779 (*. division of monomials .*)

   780 fun mv_mdiv((0,[]):mv_monom,_:mv_monom)=(0,[]):mv_monom

   781   | mv_mdiv(_,(0,[]))= raise error ("RATIONALS_MV_MDIV_EXCEPTION Division by 0 ")

   782   | mv_mdiv(p1:mv_monom,p2:mv_monom)=

   783     let

   784 	val c=ref (#1(p2));

   785 	val pp=ref [];

   786     in

   787 	(

   788 	 if !c=0 then raise error("MV_MDIV_EXCEPTION Dividing by zero")

   789 	 else c:=(#1(p1) div #1(p2));

   790 	     if #1(p2)<>0 then

   791 		 (

   792 		  pp:=(#2(mv_mmul((1,#2(p1)),(1,(map mv_minus) (#2(p2))))));

   793 		  if mv_is_negativ(!pp) then (0,!pp)

   794 		  else (!c,!pp)

   795 		      )

   796 	     else raise error("MV_MDIV_EXCEPTION Dividing by empty Polynom")

   797 		 )

   798     end;

   800 (*. prints a polynom for (internal use only) .*)

   801 fun mv_print_poly([]:mv_poly)=print("[]\n")

   802   | mv_print_poly((x,y)::[])= print("("^BasisLibrary.Int.toString(x)^","^ints2str(y)^")\n")

   803   | mv_print_poly((x,y)::p1) = (print("("^BasisLibrary.Int.toString(x)^","^ints2str(y)^"),");mv_print_poly(p1));

   806 (*. division of two multivariate polynomials .*)

   807 fun mv_division([]:mv_poly,g:mv_poly,order)=([]:mv_poly,[]:mv_poly)

   808   | mv_division(f,[],order)= raise error ("RATIONALS_MV_DIVISION_EXCEPTION Division by zero")

   809   | mv_division(f,g,order)=

   810     let

   811 	val r=ref [];

   812 	val q=ref [];

   813 	val g'=ref [];

   814 	val k=ref 0;

   815 	val m=ref (0,[0]);

   816 	val exit=ref 0;

   817     in

   818 	r := rev(sort (mv_geq order) (mv_shorten(f,order)));

   819 	g':= rev(sort (mv_geq order) (mv_shorten(g,order)));

   820 	if #1(hd(!g'))=0 then raise error("RATIONALS_MV_DIVISION_EXCEPTION: Dividing by zero") else ();

   821 	if  (mv_monom_greater (hd(!g'),hd(!r),order)) then ([(0,mv_null2(#2(hd(f))))],(!r))

   822 	else

   823 	    (

   824 	     exit:=0;

   825 	     while (if (!exit)=0 then not(mv_monom_greater (hd(!g'),hd(!r),order)) else false) do

   826 		 (

   827 		  if (#1(mv_lm(!g',order)))<>0 then m:=mv_mdiv(mv_lm(!r,order),mv_lm(!g',order))

   828 		  else raise error ("RATIONALS_MV_DIVISION_EXCEPTION: Dividing by zero");

   829 		  if #1(!m)<>0 then

   830 		      (

   831 		       q:=(!m)::(!q);

   832 		       r:=mv_add((!r),mv_minus2(mv_mul(!g',[!m],order)),order)

   833 		       )

   834 		  else exit:=1;

   835 		  if (if length(!r)<>0 then length(!g')<>0 else false) then ()

   836 		  else (exit:=1)

   837 		  );

   838 		 (rev(!q),!r)

   839 		 )

   840     end;

   842 (*. multiplies a polynomial with an integer .*)

   843 fun mv_skalar_mul([]:mv_poly,c) = []:mv_poly

   844   | mv_skalar_mul((x,y)::p1,c) = ((x * c),y)::mv_skalar_mul(p1,c);

   846 (*. inserts the a first variable into an polynomial with exponent v .*)

   847 fun mv_correct([]:mv_poly,v:int)=[]:mv_poly

   848   | mv_correct((x,y)::list,v:int)=(x,v::y)::mv_correct(list,v);

   850 (*. multivariate case .*)

   852 (*. decides if x is a factor of y .*)

   853 fun mv_divides([]:mv_poly,[]:mv_poly)=  raise error("RATIONALS_MV_DIVIDES_EXCEPTION: division by zero")

   854   | mv_divides(x,[]) =  raise error("RATIONALS_MV_DIVIDES_EXCEPTION: division by zero")

   855   | mv_divides(x:mv_poly,y:mv_poly) = #2(mv_division(y,x,LEX_))=[];

   857 (*. gets the maximum of a and b .*)

   858 fun mv_max(a,b) = if a>b then a else b;

   860 (*. gets the maximum exponent of a mv polynomial in the lexicographic term order .*)

   861 fun mv_deg([]:mv_poly) = 0

   862   | mv_deg(p1)=

   863     let

   864 	val p1'=mv_shorten(p1,LEX_);

   865     in

   866 	if length(p1')=0 then 0

   867 	else mv_max(hd(#2(hd(p1'))),mv_deg(tl(p1')))

   868     end;

   870 (*. gets the maximum exponent of a mv polynomial in the reverse lexicographic term order .*)

   871 fun mv_deg2([]:mv_poly) = 0

   872   | mv_deg2(p1)=

   873     let

   874 	val p1'=mv_shorten(p1,LEX_);

   875     in

   876 	if length(p1')=0 then 0

   877 	else mv_max(hd(rev(#2(hd(p1')))),mv_deg2(tl(p1')))

   878     end;

   880 (*. evaluates the mv polynomial at the value v of the main variable .*)

   881 fun mv_subs([]:mv_poly,v) = []:mv_poly

   882   | mv_subs((c,e)::p1:mv_poly,v) = mv_skalar_mul(mv_cut([(c,e)]),power v (hd(e))) @ mv_subs(p1,v);

   884 (*. calculates the content of a uv-polynomial in mv-representation .*)

   885 fun uv_content2([]:mv_poly) = 0

   886   | uv_content2((c,e)::p1) = (gcd_int c (uv_content2(p1)));

   888 (*. converts a uv-polynomial from mv-representation to  uv-representation .*)

   889 fun uv_to_list ([]:mv_poly)=[]:uv_poly

   890   | uv_to_list ((c1,e1)::others) =

   891     let

   892 	val count=ref 0;

   893 	val max=mv_grad((c1,e1)::others);

   894 	val help=ref ((c1,e1)::others);

   895 	val list=ref [];

   896     in

   897 	if length(e1)>1 then raise error ("RATIONALS_TO_LIST_EXCEPTION: not univariate")

   898 	else if length(e1)=0 then [c1]

   899 	     else

   900 		 (

   901 		  count:=0;

   902 		  while (!count)<=max do

   903 		      (

   904 		       if length(!help)>0 andalso hd(#2(hd(!help)))=max-(!count) then

   905 			   (

   906 			    list:=(#1(hd(!help)))::(!list);

   907 			    help:=tl(!help)

   908 			    )

   909 		       else

   910 			   (

   911 			    list:= 0::(!list)

   912 			    );

   913 		       count := (!count) + 1

   914 		       );

   915 		      (!list)

   916 		      )

   917     end;

   919 (*. converts a uv-polynomial from uv-representation to mv-representation .*)

   920 fun uv_to_poly ([]:uv_poly) = []:mv_poly

   921   | uv_to_poly p1 =

   922     let

   923 	val count=ref 0;

   924 	val help=ref p1;

   925 	val list=ref [];

   926     in

   927 	while length(!help)>0 do

   928 	    (

   929 	     if hd(!help)=0 then ()

   930 	     else list:=(hd(!help),[!count])::(!list);

   931 	     count:=(!count)+1;

   932 	     help:=tl(!help)

   933 	     );

   934 	    (!list)

   935     end;

   937 (*. univariate gcd calculation from polynomials in multivariate representation .*)

   938 fun uv_gcd ([]:mv_poly) p2 = p2

   939   | uv_gcd p1 ([]:mv_poly) = p1

   940   | uv_gcd p1 [(c,[e])] =

   941     let

   942 	val list=ref (rev(sort (mv_geq LEX_) (mv_shorten(p1,LEX_))));

   943 	val min=uv_mod_min(e,(hd(#2(hd(rev(!list))))));

   944     in

   945 	[(gcd_int (uv_content2(p1)) c,[min])]

   946     end

   947   | uv_gcd [(c,[e])] p2 =

   948     let

   949 	val list=ref (rev(sort (mv_geq LEX_) (mv_shorten(p2,LEX_))));

   950 	val min=uv_mod_min(e,(hd(#2(hd(rev(!list))))));

   951     in

   952 	[(gcd_int (uv_content2(p2)) c,[min])]

   953     end

   954   | uv_gcd p11 p22 = uv_to_poly(uv_mod_gcd (uv_to_list(mv_shorten(p11,LEX_))) (uv_to_list(mv_shorten(p22,LEX_))));

   956 (*. help function for the newton interpolation .*)

   957 fun mv_newton_help ([]:mv_poly list,k:int) = []:mv_poly list

   958   | mv_newton_help (pl:mv_poly list,k) =

   959     let

   960 	val x=ref (rev(pl));

   961 	val t=ref [];

   962 	val y=ref [];

   963 	val n=ref 1;

   964 	val n1=ref[];

   965     in

   966 	(

   967 	 while length(!x)>1 do

   968 	     (

   969 	      if length(hd(!x))>0 then n1:=mv_null2(#2(hd(hd(!x))))

   970 	      else if length(hd(tl(!x)))>0 then n1:=mv_null2(#2(hd(hd(tl(!x)))))

   971 		   else n1:=[];

   972 	      t:= #1(mv_division(mv_add(hd(!x),mv_skalar_mul(hd(tl(!x)),~1),LEX_),[(k,!n1)],LEX_));

   973 	      y:=(!t)::(!y);

   974 	      x:=tl(!x)

   975 	      );

   976 	 (!y)

   977 	 )

   978     end;

   980 (*. help function for the newton interpolation .*)

   981 fun mv_newton_add ([]:mv_poly list) t = []:mv_poly

   982   | mv_newton_add [x:mv_poly] t = x

   983   | mv_newton_add (pl:mv_poly list) t =

   984     let

   985 	val expos=ref [];

   986 	val pll=ref pl;

   987     in

   988 	(

   990 	 while length(!pll)>0 andalso hd(!pll)=[]  do

   991 	     (

   992 	      pll:=tl(!pll)

   993 	      );

   994 	 if length(!pll)>0 then expos:= #2(hd(hd(!pll))) else expos:=[];

   995 	 mv_add(hd(pl),

   996 		mv_mul(

   997 		       mv_add(mv_correct(mv_cut([(1,mv_null2(!expos))]),1),[(~t,mv_null2(!expos))],LEX_),

   998 		       mv_newton_add (tl(pl)) (t+1),

   999 		       LEX_

  1000 		       ),

  1001 		LEX_)

  1002 	 )

  1003     end;

  1005 (*. calculates the newton interpolation with polynomial coefficients .*)

  1006 (*. step-depth is 1 and if the result is not an integerpolynomial .*)

  1007 (*. this function returns [] .*)

  1008 fun mv_newton ([]:(mv_poly) list) = []:mv_poly

  1009   | mv_newton ([mp]:(mv_poly) list) = mp:mv_poly

  1010   | mv_newton pl =

  1011     let

  1012 	val c=ref pl;

  1013 	val c1=ref [];

  1014 	val n=length(pl);

  1015 	val k=ref 1;

  1016 	val i=ref n;

  1017 	val ppl=ref [];

  1018     in

  1019 	c1:=hd(pl)::[];

  1020 	c:=mv_newton_help(!c,!k);

  1021 	c1:=(hd(!c))::(!c1);

  1022 	while(length(!c)>1 andalso !k<n) do

  1023 	    (

  1024 	     k:=(!k)+1;

  1025 	     while  length(!c)>0 andalso hd(!c)=[] do c:=tl(!c);

  1026 	     if !c=[] then () else c:=mv_newton_help(!c,!k);

  1027 	     ppl:= !c;

  1028 	     if !c=[] then () else  c1:=(hd(!c))::(!c1)

  1029 	     );

  1030 	while  hd(!c1)=[] do c1:=tl(!c1);

  1031 	c1:=rev(!c1);

  1032 	ppl:= !c1;

  1033 	mv_newton_add (!c1) 1

  1034     end;

  1036 (*. sets the exponents of the first variable to zero .*)

  1037 fun mv_null3([]:mv_poly)    = []:mv_poly

  1038   | mv_null3((x,y)::xs) = (x,0::tl(y))::mv_null3(xs);

  1040 (*. calculates the minimum exponents of a multivariate polynomial .*)

  1041 fun mv_min_pp([]:mv_poly)=[]

  1042   | mv_min_pp((c,e)::xs)=

  1043     let

  1044 	val y=ref xs;

  1045 	val x=ref [];

  1046     in

  1047 	(

  1048 	 x:=e;

  1049 	 while length(!y)>0 do

  1050 	     (

  1051 	      x:=(map uv_mod_min) ((!x) ~~ (#2(hd(!y))));

  1052 	      y:=tl(!y)

  1053 	      );

  1054 	 !x

  1055 	 )

  1056     end;

  1058 (*. checks if all elements of the list have value zero .*)

  1059 fun list_is_null [] = true

  1060   | list_is_null (x::xs) = (x=0 andalso list_is_null(xs));

  1062 (* check if main variable is zero*)

  1063 fun main_zero (ms : mv_poly) = (list_is_null o (map (hd o #2))) ms;

  1065 (*. calculates the content of an polynomial .*)

  1066 fun mv_content([]:mv_poly) = []:mv_poly

  1067   | mv_content(p1) =

  1068     let

  1069 	val list=ref (rev(sort (mv_geq LEX_) (mv_shorten(p1,LEX_))));

  1070 	val test=ref (hd(#2(hd(!list))));

  1071 	val result=ref [];

  1072 	val min=(hd(#2(hd(rev(!list)))));

  1073     in

  1074 	(

  1075 	 if length(!list)>1 then

  1076 	     (

  1077 	      while (if length(!list)>0 then (hd(#2(hd(!list)))=(!test)) else false) do

  1078 		  (

  1079 		   result:=(#1(hd(!list)),tl(#2(hd(!list))))::(!result);

  1081 		   if length(!list)<1 then list:=[]

  1082 		   else list:=tl(!list)

  1084 		       );

  1085 		  if length(!list)>0 then

  1086 		   (

  1087 		    list:=mv_gcd (!result) (mv_cut(mv_content(!list)))

  1088 		    )

  1089 		  else list:=(!result);

  1090 		  list:=mv_correct(!list,0);

  1091 		  (!list)

  1092 		  )

  1093 	 else

  1094 	     (

  1095 	      mv_null3(!list)

  1096 	      )

  1097 	     )

  1098     end

  1100 (*. calculates the primitiv part of a polynomial .*)

  1101 and mv_pp([]:mv_poly) = []:mv_poly

  1102   | mv_pp(p1) = let

  1103 		    val cont=ref [];

  1104 		    val pp=ref[];

  1105 		in

  1106 		    cont:=mv_content(p1);

  1107 		    pp:=(#1(mv_division(p1,!cont,LEX_)));

  1108 		    if !pp=[]

  1109 			then raise error("RATIONALS_MV_PP_EXCEPTION: Invalid Content ")

  1110 		    else (!pp)

  1111 		end

  1113 (*. calculates the gcd of two multivariate polynomials with a modular approach .*)

  1114 and mv_gcd ([]:mv_poly) ([]:mv_poly) :mv_poly= []:mv_poly

  1115   | mv_gcd ([]:mv_poly) (p2) :mv_poly= p2:mv_poly

  1116   | mv_gcd (p1:mv_poly) ([]) :mv_poly= p1:mv_poly

  1117   | mv_gcd ([(x,xs)]:mv_poly) ([(y,ys)]):mv_poly =

  1118      let

  1119       val xpoly:mv_poly = [(x,xs)];

  1120       val ypoly:mv_poly = [(y,ys)];

  1121      in

  1122 	(

  1123 	 if xs=ys then [((gcd_int x y),xs)]

  1124 	 else [((gcd_int x y),(map uv_mod_min)(xs~~ys))]:mv_poly

  1125         )

  1126     end

  1127   | mv_gcd (p1:mv_poly) ([(y,ys)]) :mv_poly=

  1128 	(

  1129 	 [(gcd_int (uv_content2(p1)) (y),(map  uv_mod_min)(mv_min_pp(p1)~~ys))]:mv_poly

  1130 	)

  1131   | mv_gcd ([(y,ys)]:mv_poly) (p2):mv_poly =

  1132 	(

  1133          [(gcd_int (uv_content2(p2)) (y),(map  uv_mod_min)(mv_min_pp(p2)~~ys))]:mv_poly

  1134         )

  1135   | mv_gcd (p1':mv_poly) (p2':mv_poly):mv_poly=

  1136     let

  1137 	val vc=length(#2(hd(p1')));

  1138 	val cont =

  1139 		  (

  1140                    if main_zero(mv_content(p1')) andalso

  1141                      (main_zero(mv_content(p2'))) then

  1142                      mv_correct((mv_gcd (mv_cut(mv_content(p1'))) (mv_cut(mv_content(p2')))),0)

  1143                    else

  1144                      mv_gcd (mv_content(p1')) (mv_content(p2'))

  1145                   );

  1146 	val p1= #1(mv_division(p1',cont,LEX_));

  1147 	val p2= #1(mv_division(p2',cont,LEX_));

  1148 	val gcd=ref [];

  1149 	val candidate=ref [];

  1150 	val interpolation_list=ref [];

  1151 	val delta=ref [];

  1152         val p1r = ref [];

  1153         val p2r = ref [];

  1154         val p1r' = ref [];

  1155         val p2r' = ref [];

  1156 	val factor=ref [];

  1157 	val r=ref 0;

  1158 	val gcd_r=ref [];

  1159 	val d=ref 0;

  1160 	val exit=ref 0;

  1161 	val current_degree=ref 99999; (*. FIXME: unlimited ! .*)

  1162     in

  1163 	(

  1164 	 if vc<2 then (* areUnivariate(p1',p2') *)

  1165 	     (

  1166 	      gcd:=uv_gcd (mv_shorten(p1',LEX_)) (mv_shorten(p2',LEX_))

  1167 	      )

  1168 	 else

  1169 	     (

  1170 	      while !exit=0 do

  1171 		  (

  1172 		   r:=(!r)+1;

  1173                    p1r := mv_lc(p1,LEX_);

  1174 		   p2r := mv_lc(p2,LEX_);

  1175                    if main_zero(!p1r) andalso

  1176                       main_zero(!p2r)

  1177                    then

  1178                        (

  1179                         delta := mv_correct((mv_gcd (mv_cut (!p1r)) (mv_cut (!p2r))),0)

  1180                        )

  1181                    else

  1182                        (

  1183 		        delta := mv_gcd (!p1r) (!p2r)

  1184                        );

  1185 		   if mv_shorten(mv_subs(!p1r,!r),LEX_)=[] andalso

  1186 		      mv_shorten(mv_subs(!p2r,!r),LEX_)=[]

  1187                    then

  1188                        (

  1189 		       )

  1190 		   else

  1191 		       (

  1192 			gcd_r:=mv_shorten(mv_gcd (mv_shorten(mv_subs(p1,!r),LEX_))

  1193 					         (mv_shorten(mv_subs(p2,!r),LEX_)) ,LEX_);

  1194 			gcd_r:= #1(mv_division(mv_mul(mv_correct(mv_subs(!delta,!r),0),!gcd_r,LEX_),

  1195 					       mv_correct(mv_lc(!gcd_r,LEX_),0),LEX_));

  1196 			d:=mv_deg2(!gcd_r); (* deg(gcd_r,z) *)

  1197 			if (!d < !current_degree) then

  1198 			    (

  1199 			     current_degree:= !d;

  1200 			     interpolation_list:=[mv_correct(!gcd_r,0)]

  1201 			     )

  1202 			else

  1203 			    (

  1204 			     if (!d = !current_degree) then

  1205 				 (

  1206 				  interpolation_list:=mv_correct(!gcd_r,0)::(!interpolation_list)

  1207 				  )

  1208 			     else ()

  1209 				 )

  1210 			    );

  1211 		      if length(!interpolation_list)> uv_mod_min(mv_deg(p1),mv_deg(p2)) then

  1212 			  (

  1213 			   candidate := mv_newton(rev(!interpolation_list));

  1214 			   if !candidate=[] then ()

  1215 			   else

  1216 			       (

  1217 				candidate:=mv_pp(!candidate);

  1218 				if mv_divides(!candidate,p1) andalso mv_divides(!candidate,p2) then

  1219 				    (

  1220 				     gcd:= mv_mul(!candidate,cont,LEX_);

  1221 				     exit:=1

  1222 				     )

  1223 				else ()

  1224 				    );

  1225 			       interpolation_list:=[mv_correct(!gcd_r,0)]

  1226 			       )

  1227 		      else ()

  1228 			  )

  1229 	     );

  1230 	     (!gcd):mv_poly

  1231 	     )

  1232     end;

  1235 (*. calculates the least common divisor of two polynomials .*)

  1236 fun mv_lcm (p1:mv_poly) (p2:mv_poly) :mv_poly =

  1237     (

  1238      #1(mv_division(mv_mul(p1,p2,LEX_),mv_gcd p1 p2,LEX_))

  1239      );

  1241 (*. gets the variables (strings) of a term .*)

  1242 fun get_vars(term1) = (map free2str) (vars term1); (*["a","b","c"]; *)

  1244 (*. counts the negative coefficents in a polynomial .*)

  1245 fun count_neg ([]:mv_poly) = 0

  1246   | count_neg ((c,e)::xs) = if c<0 then 1+count_neg xs

  1247 			  else count_neg xs;

  1249 (*. help function for is_polynomial

  1250     checks the order of the operators .*)

  1251 fun test_polynomial (Const ("uminus",_) $ Free (str,_)) _ = true (*WN.13.3.03*)

  1252   | test_polynomial (t as Free(str,_)) v = true

  1253   | test_polynomial (t as Const ("op *",_) $ t1 $ t2) v = if v="^" then false

  1254 						     else (test_polynomial t1 "*") andalso (test_polynomial t2 "*")

  1255   | test_polynomial (t as Const ("op +",_) $ t1 $ t2) v = if v="*" orelse v="^" then false

  1256 							  else (test_polynomial t1 " ") andalso (test_polynomial t2 " ")

  1257   | test_polynomial (t as Const ("Atools.pow",_) $ t1 $ t2) v = (test_polynomial t1 "^") andalso (test_polynomial t2 "^")

  1258   | test_polynomial _ v = false;

  1260 (*. tests if a term is a polynomial .*)

  1261 fun is_polynomial t = test_polynomial t " ";

  1263 (*. help function for is_expanded

  1264     checks the order of the operators .*)

  1265 fun test_exp (t as Free(str,_)) v = true

  1266   | test_exp (t as Const ("op *",_) $ t1 $ t2) v = if v="^" then false

  1267 						     else (test_exp t1 "*") andalso (test_exp t2 "*")

  1268   | test_exp (t as Const ("op +",_) $ t1 $ t2) v = if v="*" orelse v="^" then false

  1269 							  else (test_exp t1 " ") andalso (test_exp t2 " ")

  1270   | test_exp (t as Const ("op -",_) $ t1 $ t2) v = if v="*" orelse v="^" then false

  1271 							  else (test_exp t1 " ") andalso (test_exp t2 " ")

  1272   | test_exp (t as Const ("Atools.pow",_) $ t1 $ t2) v = (test_exp t1 "^") andalso (test_exp t2 "^")

  1273   | test_exp  _ v = false;

  1276 (*. help function for check_coeff:

  1277     converts the term to a list of coefficients .*)

  1278 fun term2coef' (t as Free(str,_(*typ*))) v :mv_poly option =

  1279     let

  1280 	val x=ref None;

  1281 	val len=ref 0;

  1282 	val vl=ref [];

  1283 	val vh=ref [];

  1284 	val i=ref 0;

  1285     in

  1286 	if is_numeral str then

  1287 	    (

  1288 	     Some [(((the o int_of_str) str),mv_null2(v))] handle _ => None

  1289 		 )

  1290 	else (* variable *)

  1291 	    (

  1292 	     len:=length(v);

  1293 	     vh:=v;

  1294 	     while ((!len)>(!i)) do

  1295 		 (

  1296 		  if str=hd((!vh)) then

  1297 		      (

  1298 		       vl:=1::(!vl)

  1299 		       )

  1300 		  else

  1301 		      (

  1302 		       vl:=0::(!vl)

  1303 		       );

  1304 		      vh:=tl(!vh);

  1305 		      i:=(!i)+1

  1306 		      );

  1307 		 Some [(1,rev(!vl))] handle _ => None

  1308 	    )

  1309     end

  1310   | term2coef' (Const ("op *",_) $ t1 $ t2) v :mv_poly option=

  1311     let

  1312 	val t1pp=ref [];

  1313 	val t2pp=ref [];

  1314 	val t1c=ref 0;

  1315 	val t2c=ref 0;

  1316     in

  1317 	(

  1318 	 t1pp:=(#2(hd(the(term2coef' t1 v))));

  1319 	 t2pp:=(#2(hd(the(term2coef' t2 v))));

  1320 	 t1c:=(#1(hd(the(term2coef' t1 v))));

  1321 	 t2c:=(#1(hd(the(term2coef' t2 v))));

  1323 	 Some [( (!t1c)*(!t2c) ,( (map op+) ((!t1pp)~~(!t2pp)) ) )] handle _ => None

  1325 	 )

  1326     end

  1327   | term2coef' (Const ("Atools.pow",_) $ (t1 as Free (str1,_)) $ (t2 as Free (str2,_))) v :mv_poly option=

  1328     let

  1329 	val x=ref None;

  1330 	val len=ref 0;

  1331 	val vl=ref [];

  1332 	val vh=ref [];

  1333 	val vtemp=ref [];

  1334 	val i=ref 0;

  1335     in

  1336     (

  1337      if (not o is_numeral) str1 andalso is_numeral str2 then

  1338 	 (

  1339 	  len:=length(v);

  1340 	  vh:=v;

  1342 	  while ((!len)>(!i)) do

  1343 	      (

  1344 	       if str1=hd((!vh)) then

  1345 		   (

  1346 		    vl:=((the o int_of_str) str2)::(!vl)

  1347 		    )

  1348 	       else

  1349 		   (

  1350 		    vl:=0::(!vl)

  1351 		    );

  1352 		   vh:=tl(!vh);

  1353 		   i:=(!i)+1

  1354 		   );

  1355 	      Some [(1,rev(!vl))] handle _ => None

  1356 	      )

  1357      else raise error ("RATIONALS_TERM2COEF_EXCEPTION 1: Invalid term")

  1358 	 )

  1359     end

  1360   | term2coef' (Const ("op +",_) $ t1 $ t2) v :mv_poly option=

  1361     (

  1362      Some ((the(term2coef' t1 v)) @ (the(term2coef' t2 v))) handle _ => None

  1363 	 )

  1364   | term2coef' (Const ("op -",_) $ t1 $ t2) v :mv_poly option=

  1365     (

  1366      Some ((the(term2coef' t1 v)) @ mv_skalar_mul((the(term2coef' t2 v)),1)) handle _ => None

  1367 	 )

  1368   | term2coef' (term) v = raise error ("RATIONALS_TERM2COEF_EXCEPTION 2: Invalid term");

  1370 (*. checks if all coefficients of a polynomial are positiv (except the first) .*)

  1371 fun check_coeff t = (* erste Koeffizient kann <0 sein !!! *)

  1372     if count_neg(tl(the(term2coef' t (get_vars(t)))))=0 then true

  1373     else false;

  1375 (*. checks for expanded term [3] .*)

  1376 fun is_expanded t = test_exp t " " andalso check_coeff(t);

  1378 (*WN.7.3.03 Hilfsfunktion f"ur term2poly'*)

  1379 fun mk_monom v' p vs =

  1380     let fun conv p (v: string) = if v'= v then p else 0

  1381     in map (conv p) vs end;

  1382 (* mk_monom "y" 5 ["a","b","x","y","z"];

  1383 val it = [0,0,0,5,0] : int list*)

  1385 (*. this function converts the term representation into the internal representation mv_poly .*)

  1386 fun term2poly' (Const ("uminus",_) $ Free (str,_)) v = (*WN.7.3.03*)

  1387     if is_numeral str

  1388     then Some [((the o int_of_str) ("-"^str), mk_monom "#" 0 v)]

  1389     else Some [(~1, mk_monom str 1 v)]

  1391   | term2poly' (Free(str,_)) v :mv_poly option =

  1392     let

  1393 	val x=ref None;

  1394 	val len=ref 0;

  1395 	val vl=ref [];

  1396 	val vh=ref [];

  1397 	val i=ref 0;

  1398     in

  1399 	if is_numeral str then

  1400 	    (

  1401 	     Some [(((the o int_of_str) str),mv_null2 v)] handle _ => None

  1402 		 )

  1403 	else (* variable *)

  1404 	    (

  1405 	     len:=length v;

  1406 	     vh:= v;

  1407 	     while ((!len)>(!i)) do

  1408 		 (

  1409 		  if str=hd((!vh)) then

  1410 		      (

  1411 		       vl:=1::(!vl)

  1412 		       )

  1413 		  else

  1414 		      (

  1415 		       vl:=0::(!vl)

  1416 		       );

  1417 		      vh:=tl(!vh);

  1418 		      i:=(!i)+1

  1419 		      );

  1420 		 Some [(1,rev(!vl))] handle _ => None

  1421 	    )

  1422     end

  1423   | term2poly' (Const ("op *",_) $ t1 $ t2) v :mv_poly option=

  1424     let

  1425 	val t1pp=ref [];

  1426 	val t2pp=ref [];

  1427 	val t1c=ref 0;

  1428 	val t2c=ref 0;

  1429     in

  1430 	(

  1431 	 t1pp:=(#2(hd(the(term2poly' t1 v))));

  1432 	 t2pp:=(#2(hd(the(term2poly' t2 v))));

  1433 	 t1c:=(#1(hd(the(term2poly' t1 v))));

  1434 	 t2c:=(#1(hd(the(term2poly' t2 v))));

  1436 	 Some [( (!t1c)*(!t2c) ,( (map op+) ((!t1pp)~~(!t2pp)) ) )]

  1437 	 handle _ => None

  1439 	 )

  1440     end

  1441   | term2poly' (Const ("Atools.pow",_) $ (t1 as Free (str1,_)) $

  1442 		      (t2 as Free (str2,_))) v :mv_poly option=

  1443     let

  1444 	val x=ref None;

  1445 	val len=ref 0;

  1446 	val vl=ref [];

  1447 	val vh=ref [];

  1448 	val vtemp=ref [];

  1449 	val i=ref 0;

  1450     in

  1451     (

  1452      if (not o is_numeral) str1 andalso is_numeral str2 then

  1453 	 (

  1454 	  len:=length(v);

  1455 	  vh:=v;

  1457 	  while ((!len)>(!i)) do

  1458 	      (

  1459 	       if str1=hd((!vh)) then

  1460 		   (

  1461 		    vl:=((the o int_of_str) str2)::(!vl)

  1462 		    )

  1463 	       else

  1464 		   (

  1465 		    vl:=0::(!vl)

  1466 		    );

  1467 		   vh:=tl(!vh);

  1468 		   i:=(!i)+1

  1469 		   );

  1470 	      Some [(1,rev(!vl))] handle _ => None

  1471 	      )

  1472      else raise error ("RATIONALS_TERM2POLY_EXCEPTION 1: Invalid term")

  1473 	 )

  1474     end

  1475   | term2poly' (Const ("op +",_) $ t1 $ t2) v :mv_poly option =

  1476     (

  1477      Some ((the(term2poly' t1 v)) @ (the(term2poly' t2 v))) handle _ => None

  1478 	 )

  1479   | term2poly' (Const ("op -",_) $ t1 $ t2) v :mv_poly option =

  1480     (

  1481      Some ((the(term2poly' t1 v)) @ mv_skalar_mul((the(term2poly' t2 v)),~1)) handle _ => None

  1482 	 )

  1483   | term2poly' (term) v = raise error ("RATIONALS_TERM2POLY_EXCEPTION 2: Invalid term");

  1485 (*. translates an Isabelle term into internal representation.

  1486     term2poly

  1487     fn : term ->              (*normalform [2]                    *)

  1488     	 string list ->       (*for ...!!! BITTE DIE ERKLÄRUNG,

  1489     			       DIE DU MIR LETZTES MAL GEGEBEN HAST*)

  1490     	 mv_monom list        (*internal representation           *)

  1491     		  option      (*the translation may fail with None*)

  1492 .*)

  1493 fun term2poly (t:term) v =

  1494      if is_polynomial t then term2poly' t v

  1495      else raise error ("term2poly: invalid = "^(term2str t));

  1497 (*. same as term2poly with automatic detection of the variables .*)

  1498 fun term2polyx t = term2poly t (((map free2str) o vars) t);

  1500 (*. checks if the term is in expanded polynomial form and converts it into the internal representation .*)

  1501 fun expanded2poly (t:term) v =

  1502     (*if is_expanded t then*) term2poly' t v

  1503     (*else raise error ("RATIONALS_EXPANDED2POLY_EXCEPTION: Invalid Polynomial")*);

  1505 (*. same as expanded2poly with automatic detection of the variables .*)

  1506 fun expanded2polyx t = expanded2poly t (((map free2str) o vars) t);

  1508 (*. converts a powerproduct into term representation .*)

  1509 fun powerproduct2term(xs,v) =

  1510     let

  1511 	val xss=ref xs;

  1512 	val vv=ref v;

  1513     in

  1514 	(

  1515 	 while hd(!xss)=0 do

  1516 	     (

  1517 	      xss:=tl(!xss);

  1518 	      vv:=tl(!vv)

  1519 	      );

  1521 	 if list_is_null(tl(!xss)) then

  1522 	     (

  1523 	      if hd(!xss)=1 then Free(hd(!vv), HOLogic.realT)

  1524 	      else

  1525 		  (

  1526 		   Const("Atools.pow",[HOLogic.realT,HOLogic.realT]--->HOLogic.realT) $

  1527 		   Free(hd(!vv), HOLogic.realT) $  Free(str_of_int (hd(!xss)),HOLogic.realT)

  1528 		   )

  1529 	      )

  1530 	 else

  1531 	     (

  1532 	      if hd(!xss)=1 then

  1533 		  (

  1534 		   Const("op *",[HOLogic.realT,HOLogic.realT]--->HOLogic.realT) $

  1535 		   Free(hd(!vv), HOLogic.realT) $

  1536 		   powerproduct2term(tl(!xss),tl(!vv))

  1537 		   )

  1538 	      else

  1539 		  (

  1540 		   Const("op *",[HOLogic.realT,HOLogic.realT]--->HOLogic.realT) $

  1541 		   (

  1542 		    Const("Atools.pow",[HOLogic.realT,HOLogic.realT]--->HOLogic.realT) $

  1543 		    Free(hd(!vv), HOLogic.realT) $  Free(str_of_int (hd(!xss)),HOLogic.realT)

  1544 		    ) $

  1545 		    powerproduct2term(tl(!xss),tl(!vv))

  1546 		   )

  1547 	      )

  1548 	 )

  1549     end;

  1551 (*. converts a monom into term representation .*)

  1552 (*fun monom2term ((c,e):mv_monom, v:string list) =

  1553     if c=0 then Free(str_of_int 0,HOLogic.realT)

  1554     else

  1555 	(

  1556 	 if list_is_null(e) then

  1557 	     (

  1558 	      Free(str_of_int c,HOLogic.realT)

  1559 	      )

  1560 	 else

  1561 	     (

  1562 	      if c=1 then

  1563 		  (

  1564 		   powerproduct2term(e,v)

  1565 		   )

  1566 	      else

  1567 		  (

  1568 		   Const("op *",[HOLogic.realT,HOLogic.realT]--->HOLogic.realT) $

  1569 		   Free(str_of_int c,HOLogic.realT)  $

  1570 		   powerproduct2term(e,v)

  1571 		   )

  1572 		  )

  1573 	     );*)

  1576 (*fun monom2term ((i, is):mv_monom, v) =

  1577     if list_is_null is

  1578     then

  1579 	if i >= 0

  1580 	then Free (str_of_int i, HOLogic.realT)

  1581 	else Const ("uminus", HOLogic.realT --> HOLogic.realT) $

  1582 		   Free ((str_of_int o abs) i, HOLogic.realT)

  1583     else

  1584 	if i > 0

  1585 	then Const ("op *", [HOLogic.realT,HOLogic.realT]---> HOLogic.realT) $

  1586 		   (Free (str_of_int i, HOLogic.realT)) $

  1587 		   powerproduct2term(is, v)

  1588 	else Const ("op *", [HOLogic.realT,HOLogic.realT]---> HOLogic.realT) $

  1589 		   (Const ("uminus", HOLogic.realT --> HOLogic.realT) $

  1590 		   Free ((str_of_int o abs) i, HOLogic.realT)) $

  1591 		   powerproduct2term(is, vs);---------------------------*)

  1592 fun monom2term ((i, is) : mv_monom, vs) =

  1593     if list_is_null is

  1594     then Free (str_of_int i, HOLogic.realT)

  1595     else if i = 1

  1596     then powerproduct2term (is, vs)

  1597     else Const ("op *", [HOLogic.realT, HOLogic.realT] ---> HOLogic.realT) $

  1598 	       (Free (str_of_int i, HOLogic.realT)) $

  1599 	       powerproduct2term (is, vs);

  1601 (*. converts the internal polynomial representation into an Isabelle term.*)

  1602 fun poly2term' ([] : mv_poly, vs) = Free(str_of_int 0, HOLogic.realT)

  1603   | poly2term' ([(c, e) : mv_monom], vs) = monom2term ((c, e), vs)

  1604   | poly2term' ((c, e) :: ces, vs) =

  1605     Const("op +", [HOLogic.realT, HOLogic.realT] ---> HOLogic.realT) $

  1606          poly2term (ces, vs) $ monom2term ((c, e), vs)

  1607 and poly2term (xs, vs) = poly2term' (rev (sort (mv_geq LEX_) (xs)), vs);

  1610 (*. converts a monom into term representation .*)

  1611 (*. ignores the sign of the coefficients => use only for exp-poly functions .*)

  1612 fun monom2term2((c,e):mv_monom, v:string list) =

  1613     if c=0 then Free(str_of_int 0,HOLogic.realT)

  1614     else

  1615 	(

  1616 	 if list_is_null(e) then

  1617 	     (

  1618 	      Free(str_of_int (abs(c)),HOLogic.realT)

  1619 	      )

  1620 	 else

  1621 	     (

  1622 	      if abs(c)=1 then

  1623 		  (

  1624 		   powerproduct2term(e,v)

  1625 		   )

  1626 	      else

  1627 		  (

  1628 		   Const("op *",[HOLogic.realT,HOLogic.realT]--->HOLogic.realT) $

  1629 		   Free(str_of_int (abs(c)),HOLogic.realT)  $

  1630 		   powerproduct2term(e,v)

  1631 		   )

  1632 		  )

  1633 	     );

  1635 (*. converts the expanded polynomial representation into the term representation .*)

  1636 fun exp2term' ([]:mv_poly,vars) =  Free(str_of_int 0,HOLogic.realT)

  1637   | exp2term' ([(c,e)],vars) =     monom2term((c,e),vars)

  1638   | exp2term' ((c1,e1)::others,vars) =

  1639     if c1<0 then

  1640 	Const("op -",[HOLogic.realT,HOLogic.realT]--->HOLogic.realT) $

  1641 	exp2term'(others,vars) $

  1642 	(

  1643 	 monom2term2((c1,e1),vars)

  1644 	 )

  1645     else

  1646 	Const("op +",[HOLogic.realT,HOLogic.realT]--->HOLogic.realT) $

  1647 	exp2term'(others,vars) $

  1648 	(

  1649 	 monom2term2((c1,e1),vars)

  1650 	 );

  1652 (*. sorts the powerproduct by lexicographic termorder and converts them into

  1653     a term in polynomial representation .*)

  1654 fun poly2expanded (xs,vars) = exp2term'(rev(sort (mv_geq LEX_) (xs)),vars);

  1656 (*. converts a polynomial into expanded form .*)

  1657 fun polynomial2expanded t =

  1658     (let

  1659 	val vars=(((map free2str) o vars) t);

  1660     in

  1661 	Some (poly2expanded (the (term2poly t vars), vars))

  1662     end) handle _ => None;

  1664 (*. converts a polynomial into polynomial form .*)

  1665 fun expanded2polynomial t =

  1666     (let

  1667 	val vars=(((map free2str) o vars) t);

  1668     in

  1669 	Some (poly2term (the (expanded2poly t vars), vars))

  1670     end) handle _ => None;

  1673 (*. calculates the greatest common divisor of numerator and denominator and seperates it from each .*)

  1674 fun step_cancel (t as Const ("HOL.divide",_) $ p1 $ p2) =

  1675     let

  1676 	val p1' = ref [];

  1677 	val p2' = ref [];

  1678 	val p3  = ref []

  1679 	val vars = rev(get_vars(p1) union get_vars(p2));

  1680     in

  1681 	(

  1682          p1':= sort (mv_geq LEX_) (the (term2poly p1 vars ));

  1683        	 p2':= sort (mv_geq LEX_) (the (term2poly p2 vars ));

  1684 	 p3:= sort (mv_geq LEX_) (mv_gcd (!p1') (!p2'));

  1685 	 if (!p3)=[(1,mv_null2(vars))] then

  1686 	     (

  1687 	      Const ("HOL.divide",[HOLogic.realT,HOLogic.realT]--->HOLogic.realT) $ p1 $ p2

  1688 	      )

  1689 	 else

  1690 	     (

  1692 	      p1':=sort (mv_geq LEX_) (#1(mv_division((!p1'),(!p3),LEX_)));

  1693 	      p2':=sort (mv_geq LEX_) (#1(mv_division((!p2'),(!p3),LEX_)));

  1695 	      if #1(hd(sort (mv_geq LEX_) (!p2'))) (*mv_lc2(!p2',LEX_)*)>0 then

  1696 	      (

  1697 	       Const ("HOL.divide",[HOLogic.realT,HOLogic.realT]--->HOLogic.realT)

  1698 	       $

  1699 	       (

  1700 		Const ("op *",[HOLogic.realT,HOLogic.realT]--->HOLogic.realT) $

  1701 		poly2term(!p1',vars) $

  1702 		poly2term(!p3,vars)

  1703 		)

  1704 	       $

  1705 	       (

  1706 		Const ("op *",[HOLogic.realT,HOLogic.realT]--->HOLogic.realT) $

  1707 		poly2term(!p2',vars) $

  1708 		poly2term(!p3,vars)

  1709 		)

  1710 	       )

  1711 	      else

  1712 	      (

  1713 	       p1':=mv_skalar_mul(!p1',~1);

  1714 	       p2':=mv_skalar_mul(!p2',~1);

  1715 	       p3:=mv_skalar_mul(!p3,~1);

  1716 	       (

  1717 		Const ("HOL.divide",[HOLogic.realT,HOLogic.realT]--->HOLogic.realT)

  1718 		$

  1719 		(

  1720 		 Const ("op *",[HOLogic.realT,HOLogic.realT]--->HOLogic.realT) $

  1721 		 poly2term(!p1',vars) $

  1722 		 poly2term(!p3,vars)

  1723 		 )

  1724 		$

  1725 		(

  1726 		 Const ("op *",[HOLogic.realT,HOLogic.realT]--->HOLogic.realT) $

  1727 		 poly2term(!p2',vars) $

  1728 		 poly2term(!p3,vars)

  1729 		 )

  1730 		)

  1731 	       )

  1732 	      )

  1733 	     )

  1734     end

  1735 | step_cancel _ = raise error ("RATIONALS_STEP_CANCEL_EXCEPTION: Invalid fraction");

  1738 (*. same as step_cancel, this time for expanded forms (input+output) .*)

  1739 fun step_cancel_expanded (t as Const ("HOL.divide",_) $ p1 $ p2) =

  1740     let

  1741 	val p1' = ref [];

  1742 	val p2' = ref [];

  1743 	val p3  = ref []

  1744 	val vars = rev(get_vars(p1) union get_vars(p2));

  1745     in

  1746 	(

  1747          p1':= sort (mv_geq LEX_) (the (expanded2poly p1 vars ));

  1748        	 p2':= sort (mv_geq LEX_) (the (expanded2poly p2 vars ));

  1749 	 p3:= sort (mv_geq LEX_) (mv_gcd (!p1') (!p2'));

  1750 	 if (!p3)=[(1,mv_null2(vars))] then

  1751 	     (

  1752 	      Const ("HOL.divide",[HOLogic.realT,HOLogic.realT]--->HOLogic.realT) $ p1 $ p2

  1753 	      )

  1754 	 else

  1755 	     (

  1757 	      p1':=sort (mv_geq LEX_) (#1(mv_division((!p1'),(!p3),LEX_)));

  1758 	      p2':=sort (mv_geq LEX_) (#1(mv_division((!p2'),(!p3),LEX_)));

  1760 	      if #1(hd(sort (mv_geq LEX_) (!p2')))(* mv_lc2(!p2',LEX_)*)>0 then

  1761 	      (

  1762 	       Const ("HOL.divide",[HOLogic.realT,HOLogic.realT]--->HOLogic.realT)

  1763 	       $

  1764 	       (

  1765 		Const ("op *",[HOLogic.realT,HOLogic.realT]--->HOLogic.realT) $

  1766 		poly2expanded(!p1',vars) $

  1767 		poly2expanded(!p3,vars)

  1768 		)

  1769 	       $

  1770 	       (

  1771 		Const ("op *",[HOLogic.realT,HOLogic.realT]--->HOLogic.realT) $

  1772 		poly2expanded(!p2',vars) $

  1773 		poly2expanded(!p3,vars)

  1774 		)

  1775 	       )

  1776 	      else

  1777 	      (

  1778 	       p1':=mv_skalar_mul(!p1',~1);

  1779 	       p2':=mv_skalar_mul(!p2',~1);

  1780 	       p3:=mv_skalar_mul(!p3,~1);

  1781 	       (

  1782 		Const ("HOL.divide",[HOLogic.realT,HOLogic.realT]--->HOLogic.realT)

  1783 		$

  1784 		(

  1785 		 Const ("op *",[HOLogic.realT,HOLogic.realT]--->HOLogic.realT) $

  1786 		 poly2expanded(!p1',vars) $

  1787 		 poly2expanded(!p3,vars)

  1788 		 )

  1789 		$

  1790 		(

  1791 		 Const ("op *",[HOLogic.realT,HOLogic.realT]--->HOLogic.realT) $

  1792 		 poly2expanded(!p2',vars) $

  1793 		 poly2expanded(!p3,vars)

  1794 		 )

  1795 		)

  1796 	       )

  1797 	      )

  1798 	     )

  1799     end

  1800 | step_cancel_expanded _ = raise error ("RATIONALS_STEP_CANCEL_EXCEPTION: Invalid fraction");

  1802 (*. calculates the greatest common divisor of numerator and denominator and divides each through it .*)

  1803 fun direct_cancel (t as Const ("HOL.divide",_) $ p1 $ p2) =

  1804     let

  1805 	val p1' = ref [];

  1806 	val p2' = ref [];

  1807 	val p3  = ref []

  1808 	val vars = rev(get_vars(p1) union get_vars(p2));

  1809     in

  1810 	(

  1811 	 p1':=sort (mv_geq LEX_) (mv_shorten((the (term2poly p1 vars )),LEX_));

  1812 	 p2':=sort (mv_geq LEX_) (mv_shorten((the (term2poly p2 vars )),LEX_));

  1813 	 p3 :=sort (mv_geq LEX_) (mv_gcd (!p1') (!p2'));

  1815 	 if (!p3)=[(1,mv_null2(vars))] then

  1816 	     (

  1817 	      (Const ("HOL.divide",[HOLogic.realT,HOLogic.realT]--->HOLogic.realT) $ p1 $ p2,[])

  1818 	      )

  1819 	 else

  1820 	     (

  1821 	      p1':=sort (mv_geq LEX_) (#1(mv_division((!p1'),(!p3),LEX_)));

  1822 	      p2':=sort (mv_geq LEX_) (#1(mv_division((!p2'),(!p3),LEX_)));

  1823 	      if #1(hd(sort (mv_geq LEX_) (!p2'))) (*mv_lc2(!p2',LEX_)*)>0 then

  1824 	      (

  1825 	       (

  1826 		Const ("HOL.divide",[HOLogic.realT,HOLogic.realT]--->HOLogic.realT)

  1827 		$

  1828 		(

  1829 		 poly2term((!p1'),vars)

  1830 		 )

  1831 		$

  1832 		(

  1833 		 poly2term((!p2'),vars)

  1834 		 )

  1835 		)

  1836 	       ,

  1837 	       if mv_grad(!p3)>0 then

  1838 		   [

  1839 		    (

  1840 		     Const ("Not",[bool]--->bool) $

  1841 		     (

  1842 		      Const("op =",[HOLogic.realT,HOLogic.realT]--->bool) $

  1843 		      poly2term((!p3),vars) $

  1844 		      Free("0",HOLogic.realT)

  1845 		      )

  1846 		     )

  1847 		    ]

  1848 	       else

  1849 		   []

  1850 		   )

  1851 	      else

  1852 		  (

  1853 		   p1':=mv_skalar_mul(!p1',~1);

  1854 		   p2':=mv_skalar_mul(!p2',~1);

  1855 		   if length(!p3)> 2*(count_neg(!p3)) then () else p3 :=mv_skalar_mul(!p3,~1);

  1856 		       (

  1857 			Const ("HOL.divide",[HOLogic.realT,HOLogic.realT]--->HOLogic.realT)

  1858 			$

  1859 			(

  1860 			 poly2term((!p1'),vars)

  1861 			 )

  1862 			$

  1863 			(

  1864 			 poly2term((!p2'),vars)

  1865 			 )

  1866 			,

  1867 			if mv_grad(!p3)>0 then

  1868 			    [

  1869 			     (

  1870 			      Const ("Not",[bool]--->bool) $

  1871 			      (

  1872 			       Const("op =",[HOLogic.realT,HOLogic.realT]--->bool) $

  1873 			       poly2term((!p3),vars) $

  1874 			       Free("0",HOLogic.realT)

  1875 			       )

  1876 			      )

  1877 			     ]

  1878 			else

  1879 			    []

  1880 			    )

  1881 		       )

  1882 		  )

  1883 	     )

  1884     end

  1885   | direct_cancel _ = raise error ("RATIONALS_DIRECT_CANCEL_EXCEPTION: Invalid fraction");

  1887 (*. same es direct_cancel, this time for expanded forms (input+output).*)

  1888 fun direct_cancel_expanded (t as Const ("HOL.divide",_) $ p1 $ p2) =

  1889     let

  1890 	val p1' = ref [];

  1891 	val p2' = ref [];

  1892 	val p3  = ref []

  1893 	val vars = rev(get_vars(p1) union get_vars(p2));

  1894     in

  1895 	(

  1896 	 p1':=sort (mv_geq LEX_) (mv_shorten((the (expanded2poly p1 vars )),LEX_));

  1897 	 p2':=sort (mv_geq LEX_) (mv_shorten((the (expanded2poly p2 vars )),LEX_));

  1898 	 p3 :=sort (mv_geq LEX_) (mv_gcd (!p1') (!p2'));

  1900 	 if (!p3)=[(1,mv_null2(vars))] then

  1901 	     (

  1902 	      (Const ("HOL.divide",[HOLogic.realT,HOLogic.realT]--->HOLogic.realT) $ p1 $ p2,[])

  1903 	      )

  1904 	 else

  1905 	     (

  1906 	      p1':=sort (mv_geq LEX_) (#1(mv_division((!p1'),(!p3),LEX_)));

  1907 	      p2':=sort (mv_geq LEX_) (#1(mv_division((!p2'),(!p3),LEX_)));

  1908 	      if #1(hd(sort (mv_geq LEX_) (!p2'))) (*mv_lc2(!p2',LEX_)*)>0 then

  1909 	      (

  1910 	       (

  1911 		Const ("HOL.divide",[HOLogic.realT,HOLogic.realT]--->HOLogic.realT)

  1912 		$

  1913 		(

  1914 		 poly2expanded((!p1'),vars)

  1915 		 )

  1916 		$

  1917 		(

  1918 		 poly2expanded((!p2'),vars)

  1919 		 )

  1920 		)

  1921 	       ,

  1922 	       if mv_grad(!p3)>0 then

  1923 		   [

  1924 		    (

  1925 		     Const ("Not",[bool]--->bool) $

  1926 		     (

  1927 		      Const("op =",[HOLogic.realT,HOLogic.realT]--->bool) $

  1928 		      poly2expanded((!p3),vars) $

  1929 		      Free("0",HOLogic.realT)

  1930 		      )

  1931 		     )

  1932 		    ]

  1933 	       else

  1934 		   []

  1935 		   )

  1936 	      else

  1937 		  (

  1938 		   p1':=mv_skalar_mul(!p1',~1);

  1939 		   p2':=mv_skalar_mul(!p2',~1);

  1940 		   if length(!p3)> 2*(count_neg(!p3)) then () else p3 :=mv_skalar_mul(!p3,~1);

  1941 		       (

  1942 			Const ("HOL.divide",[HOLogic.realT,HOLogic.realT]--->HOLogic.realT)

  1943 			$

  1944 			(

  1945 			 poly2expanded((!p1'),vars)

  1946 			 )

  1947 			$

  1948 			(

  1949 			 poly2expanded((!p2'),vars)

  1950 			 )

  1951 			,

  1952 			if mv_grad(!p3)>0 then

  1953 			    [

  1954 			     (

  1955 			      Const ("Not",[bool]--->bool) $

  1956 			      (

  1957 			       Const("op =",[HOLogic.realT,HOLogic.realT]--->bool) $

  1958 			       poly2expanded((!p3),vars) $

  1959 			       Free("0",HOLogic.realT)

  1960 			       )

  1961 			      )

  1962 			     ]

  1963 			else

  1964 			    []

  1965 			    )

  1966 		       )

  1967 		  )

  1968 	     )

  1969     end

  1970   | direct_cancel_expanded _ = raise error ("RATIONALS_DIRECT_CANCEL_EXCEPTION: Invalid fraction");

  1973 (*. adds two fractions .*)

  1974 fun add_fract ((Const("HOL.divide",_) $ t11 $ t12),(Const("HOL.divide",_) $ t21 $ t22)) =

  1975     let

  1976 	val vars=get_vars(t11) union get_vars(t12) union get_vars(t21) union get_vars(t22);

  1977 	val t11'=ref (the(term2poly t11 vars));

  1978 val _= writeln"### add_fract: done t11"

  1979 	val t12'=ref (the(term2poly t12 vars));

  1980 val _= writeln"### add_fract: done t12"

  1981 	val t21'=ref (the(term2poly t21 vars));

  1982 val _= writeln"### add_fract: done t21"

  1983 	val t22'=ref (the(term2poly t22 vars));

  1984 val _= writeln"### add_fract: done t22"

  1985 	val den=ref [];

  1986 	val nom=ref [];

  1987 	val m1=ref [];

  1988 	val m2=ref [];

  1989     in

  1991 	(

  1992 	 den :=sort (mv_geq LEX_) (mv_lcm (!t12') (!t22'));

  1993 writeln"### add_fract: done sort mv_lcm";

  1994 	 m1  :=sort (mv_geq LEX_) (#1(mv_division(!den,!t12',LEX_)));

  1995 writeln"### add_fract: done sort mv_division t12";

  1996 	 m2  :=sort (mv_geq LEX_) (#1(mv_division(!den,!t22',LEX_)));

  1997 writeln"### add_fract: done sort mv_division t22";

  1998 	 nom :=sort (mv_geq LEX_)

  1999 		    (mv_shorten(mv_add(mv_mul(!t11',!m1,LEX_),

  2000 				       mv_mul(!t21',!m2,LEX_),

  2001 				       LEX_),

  2002 				LEX_));

  2003 writeln"### add_fract: done sort mv_add";

  2004 	 (

  2005 	  Const ("HOL.divide",[HOLogic.realT,HOLogic.realT]--->HOLogic.realT)

  2006 	  $

  2007 	  (

  2008 	   poly2term((!nom),vars)

  2009 	   )

  2010 	  $

  2011 	  (

  2012 	   poly2term((!den),vars)

  2013 	   )

  2014 	  )

  2015 	 )

  2016     end

  2017   | add_fract (_,_) = raise error ("RATIONALS_ADD_FRACTION_EXCEPTION: Invalid add_fraction call");

  2019 (*. adds two expanded fractions .*)

  2020 fun add_fract_exp ((Const("HOL.divide",_) $ t11 $ t12),(Const("HOL.divide",_) $ t21 $ t22)) =

  2021     let

  2022 	val vars=get_vars(t11) union get_vars(t12) union get_vars(t21) union get_vars(t22);

  2023 	val t11'=ref (the(expanded2poly t11 vars));

  2024 	val t12'=ref (the(expanded2poly t12 vars));

  2025 	val t21'=ref (the(expanded2poly t21 vars));

  2026 	val t22'=ref (the(expanded2poly t22 vars));

  2027 	val den=ref [];

  2028 	val nom=ref [];

  2029 	val m1=ref [];

  2030 	val m2=ref [];

  2031     in

  2033 	(

  2034 	 den :=sort (mv_geq LEX_) (mv_lcm (!t12') (!t22'));

  2035 	 m1  :=sort (mv_geq LEX_) (#1(mv_division(!den,!t12',LEX_)));

  2036 	 m2  :=sort (mv_geq LEX_) (#1(mv_division(!den,!t22',LEX_)));

  2037 	 nom :=sort (mv_geq LEX_) (mv_shorten(mv_add(mv_mul(!t11',!m1,LEX_),mv_mul(!t21',!m2,LEX_),LEX_),LEX_));

  2038 	 (

  2039 	  Const ("HOL.divide",[HOLogic.realT,HOLogic.realT]--->HOLogic.realT)

  2040 	  $

  2041 	  (

  2042 	   poly2expanded((!nom),vars)

  2043 	   )

  2044 	  $

  2045 	  (

  2046 	   poly2expanded((!den),vars)

  2047 	   )

  2048 	  )

  2049 	 )

  2050     end

  2051   | add_fract_exp (_,_) = raise error ("RATIONALS_ADD_FRACTION_EXP_EXCEPTION: Invalid add_fraction call");

  2053 (*. adds a list of terms .*)

  2054 fun add_list_of_fractions []= (Free("0",HOLogic.realT),[])

  2055   | add_list_of_fractions [x]= direct_cancel x

  2056   | add_list_of_fractions (x::y::xs) =

  2057     let

  2058 	val (t1a,rest1)=direct_cancel(x);

  2059 val _= writeln"### add_list_of_fractions xs: has done direct_cancel(x)";

  2060 	val (t2a,rest2)=direct_cancel(y);

  2061 val _= writeln"### add_list_of_fractions xs: has done direct_cancel(y)";

  2062 	val (t3a,rest3)=(add_list_of_fractions (add_fract(t1a,t2a)::xs));

  2063 val _= writeln"### add_list_of_fractions xs: has done add_list_of_fraction xs";

  2064 	val (t4a,rest4)=direct_cancel(t3a);

  2065 val _= writeln"### add_list_of_fractions xs: has done direct_cancel(t3a)";

  2066 	val rest=rest1 union rest2 union rest3 union rest4;

  2067     in

  2068 	(writeln"### add_list_of_fractions in";

  2069 	 (

  2070 	 (t4a,rest)

  2071 	 )

  2072 	 )

  2073     end;

  2075 (*. adds a list of expanded terms .*)

  2076 fun add_list_of_fractions_exp []= (Free("0",HOLogic.realT),[])

  2077   | add_list_of_fractions_exp [x]= direct_cancel_expanded x

  2078   | add_list_of_fractions_exp (x::y::xs) =

  2079     let

  2080 	val (t1a,rest1)=direct_cancel_expanded(x);

  2081 	val (t2a,rest2)=direct_cancel_expanded(y);

  2082 	val (t3a,rest3)=(add_list_of_fractions_exp (add_fract_exp(t1a,t2a)::xs));

  2083 	val (t4a,rest4)=direct_cancel_expanded(t3a);

  2084 	val rest=rest1 union rest2 union rest3 union rest4;

  2085     in

  2086 	(

  2087 	 (t4a,rest)

  2088 	 )

  2089     end;

  2091 (*. calculates the lcm of a list of mv_poly .*)

  2092 fun calc_lcm ([x],var)= (x,var)

  2093   | calc_lcm ((x::xs),var) = (mv_lcm x (#1(calc_lcm (xs,var))),var);

  2095 (*. converts a list of terms to a list of mv_poly .*)

  2096 fun t2d([],_)=[]

  2097   | t2d((t as (Const("HOL.divide",_) $ p1 $ p2))::xs,vars)= (the(term2poly p2 vars)) :: t2d(xs,vars);

  2099 (*. same as t2d, this time for expanded forms .*)

  2100 fun t2d_exp([],_)=[]

  2101   | t2d_exp((t as (Const("HOL.divide",_) $ p1 $ p2))::xs,vars)= (the(expanded2poly p2 vars)) :: t2d_exp(xs,vars);

  2103 (*. converts a list of fract terms to a list of their denominators .*)

  2104 fun termlist2denominators [] = ([],[])

  2105   | termlist2denominators xs =

  2106     let

  2107 	val xxs=ref xs;

  2108 	val var=ref [];

  2109     in

  2110 	var:=[];

  2111 	while length(!xxs)>0 do

  2112 	    (

  2113 	     let

  2114 		 val (t as Const ("HOL.divide",_) $ p1x $ p2x)=hd(!xxs);

  2115 	     in

  2116 		 (

  2117 		  xxs:=tl(!xxs);

  2118 		  var:=((get_vars(p2x)) union (get_vars(p1x)) union (!var))

  2119 		  )

  2120 	     end

  2121 	     );

  2122 	    (t2d(xs,!var),!var)

  2123     end;

  2125 (*. calculates the lcm of a list of mv_poly .*)

  2126 fun calc_lcm ([x],var)= (x,var)

  2127   | calc_lcm ((x::xs),var) = (mv_lcm x (#1(calc_lcm (xs,var))),var);

  2129 (*. converts a list of terms to a list of mv_poly .*)

  2130 fun t2d([],_)=[]

  2131   | t2d((t as (Const("HOL.divide",_) $ p1 $ p2))::xs,vars)= (the(term2poly p2 vars)) :: t2d(xs,vars);

  2133 (*. same as t2d, this time for expanded forms .*)

  2134 fun t2d_exp([],_)=[]

  2135   | t2d_exp((t as (Const("HOL.divide",_) $ p1 $ p2))::xs,vars)= (the(expanded2poly p2 vars)) :: t2d_exp(xs,vars);

  2137 (*. converts a list of fract terms to a list of their denominators .*)

  2138 fun termlist2denominators [] = ([],[])

  2139   | termlist2denominators xs =

  2140     let

  2141 	val xxs=ref xs;

  2142 	val var=ref [];

  2143     in

  2144 	var:=[];

  2145 	while length(!xxs)>0 do

  2146 	    (

  2147 	     let

  2148 		 val (t as Const ("HOL.divide",_) $ p1x $ p2x)=hd(!xxs);

  2149 	     in

  2150 		 (

  2151 		  xxs:=tl(!xxs);

  2152 		  var:=((get_vars(p2x)) union (get_vars(p1x)) union (!var))

  2153 		  )

  2154 	     end

  2155 	     );

  2156 	    (t2d(xs,!var),!var)

  2157     end;

  2159 (*. same as termlist2denminators, this time for expanded forms .*)

  2160 fun termlist2denominators_exp [] = ([],[])

  2161   | termlist2denominators_exp xs =

  2162     let

  2163 	val xxs=ref xs;

  2164 	val var=ref [];

  2165     in

  2166 	var:=[];

  2167 	while length(!xxs)>0 do

  2168 	    (

  2169 	     let

  2170 		 val (t as Const ("HOL.divide",_) $ p1x $ p2x)=hd(!xxs);

  2171 	     in

  2172 		 (

  2173 		  xxs:=tl(!xxs);

  2174 		  var:=((get_vars(p2x)) union (get_vars(p1x)) union (!var))

  2175 		  )

  2176 	     end

  2177 	     );

  2178 	    (t2d_exp(xs,!var),!var)

  2179     end;

  2181 (*. reduces all fractions to the least common denominator .*)

  2182 fun com_den(x::xs,denom,den,var)=

  2183     let

  2184 	val (t as Const ("HOL.divide",_) $ p1' $ p2')=x;

  2185 	val p2= sort (mv_geq LEX_) (the(term2poly p2' var));

  2186 	val p3= #1(mv_division(denom,p2,LEX_));

  2187 	val p1var=get_vars(p1');

  2188     in

  2189 	if length(xs)>0 then

  2190 	    if p3=[(1,mv_null2(var))] then

  2191 		(

  2192 		 Const ("op +",[HOLogic.realT,HOLogic.realT]--->HOLogic.realT)

  2193 		 $

  2194 		 (

  2195 		  Const ("HOL.divide",[HOLogic.realT,HOLogic.realT]--->HOLogic.realT)

  2196 		  $

  2197 		  poly2term(the (term2poly p1' p1var),p1var)

  2198 		  $

  2199 		  den

  2200 		  )

  2201 		 $

  2202 		 #1(com_den(xs,denom,den,var))

  2203 		,

  2204 		[]

  2205 		)

  2206 	    else

  2207 		(

  2208 		 Const ("op +",[HOLogic.realT,HOLogic.realT]--->HOLogic.realT)

  2209 		 $

  2210 		 (

  2211 		  Const ("HOL.divide",[HOLogic.realT,HOLogic.realT]--->HOLogic.realT)

  2212 		  $

  2213 		  (

  2214 		   Const ("op *",[HOLogic.realT,HOLogic.realT]--->HOLogic.realT) $

  2215 		   poly2term(the (term2poly p1' p1var),p1var) $

  2216 		   poly2term(p3,var)

  2217 		   )

  2218 		  $

  2219 		  (

  2220 		   den

  2221 		   )

  2222 		  )

  2223 		 $

  2224 		 #1(com_den(xs,denom,den,var))

  2225 		,

  2226 		[]

  2227 		)

  2228 	else

  2229 	    if p3=[(1,mv_null2(var))] then

  2230 		(

  2231 		 (

  2232 		  Const ("HOL.divide",[HOLogic.realT,HOLogic.realT]--->HOLogic.realT)

  2233 		  $

  2234 		  poly2term(the (term2poly p1' p1var),p1var)

  2235 		  $

  2236 		  den

  2237 		  )

  2238 		 ,

  2239 		 []

  2240 		 )

  2241 	     else

  2242 		 (

  2243 		  Const ("HOL.divide",[HOLogic.realT,HOLogic.realT]--->HOLogic.realT)

  2244 		  $

  2245 		  (

  2246 		   Const ("op *",[HOLogic.realT,HOLogic.realT]--->HOLogic.realT) $

  2247 		   poly2term(the (term2poly p1' p1var),p1var) $

  2248 		   poly2term(p3,var)

  2249 		   )

  2250 		  $

  2251 		  den

  2252 		  ,

  2253 		  []

  2254 		  )

  2255     end;

  2257 (*. same as com_den, this time for expanded forms .*)

  2258 fun com_den_exp(x::xs,denom,den,var)=

  2259     let

  2260 	val (t as Const ("HOL.divide",_) $ p1' $ p2')=x;

  2261 	val p2= sort (mv_geq LEX_) (the(expanded2poly p2' var));

  2262 	val p3= #1(mv_division(denom,p2,LEX_));

  2263 	val p1var=get_vars(p1');

  2264     in

  2265 	if length(xs)>0 then

  2266 	    if p3=[(1,mv_null2(var))] then

  2267 		(

  2268 		 Const ("op +",[HOLogic.realT,HOLogic.realT]--->HOLogic.realT)

  2269 		 $

  2270 		 (

  2271 		  Const ("HOL.divide",[HOLogic.realT,HOLogic.realT]--->HOLogic.realT)

  2272 		  $

  2273 		  poly2expanded(the(expanded2poly p1' p1var),p1var)

  2274 		  $

  2275 		  den

  2276 		  )

  2277 		 $

  2278 		 #1(com_den_exp(xs,denom,den,var))

  2279 		,

  2280 		[]

  2281 		)

  2282 	    else

  2283 		(

  2284 		 Const ("op +",[HOLogic.realT,HOLogic.realT]--->HOLogic.realT)

  2285 		 $

  2286 		 (

  2287 		  Const ("HOL.divide",[HOLogic.realT,HOLogic.realT]--->HOLogic.realT)

  2288 		  $

  2289 		  (

  2290 		   Const ("op *",[HOLogic.realT,HOLogic.realT]--->HOLogic.realT) $

  2291 		   poly2expanded(the(expanded2poly p1' p1var),p1var) $

  2292 		   poly2expanded(p3,var)

  2293 		   )

  2294 		  $

  2295 		  (

  2296 		   den

  2297 		   )

  2298 		  )

  2299 		 $

  2300 		 #1(com_den_exp(xs,denom,den,var))

  2301 		,

  2302 		[]

  2303 		)

  2304 	else

  2305 	    if p3=[(1,mv_null2(var))] then

  2306 		(

  2307 		 (

  2308 		  Const ("HOL.divide",[HOLogic.realT,HOLogic.realT]--->HOLogic.realT)

  2309 		  $

  2310 		  poly2expanded(the(expanded2poly p1' p1var),p1var)

  2311 		  $

  2312 		  den

  2313 		  )

  2314 		 ,

  2315 		 []

  2316 		 )

  2317 	     else

  2318 		 (

  2319 		  Const ("HOL.divide",[HOLogic.realT,HOLogic.realT]--->HOLogic.realT)

  2320 		  $

  2321 		  (

  2322 		   Const ("op *",[HOLogic.realT,HOLogic.realT]--->HOLogic.realT) $

  2323 		   poly2expanded(the(expanded2poly p1' p1var),p1var) $

  2324 		   poly2expanded(p3,var)

  2325 		   )

  2326 		  $

  2327 		  den

  2328 		  ,

  2329 		  []

  2330 		  )

  2331     end;

  2333 (* wird aktuell nicht mehr gebraucht, bei rückänderung schon

  2334 -------------------------------------------------------------

  2335 (* WN0210???SK brauch ma des überhaupt *)

  2336 fun com_den2(x::xs,denom,den,var)=

  2337     let

  2338 	val (t as Const ("HOL.divide",_) $ p1' $ p2')=x;

  2339 	val p2= sort (mv_geq LEX_) (the(term2poly p2' var));

  2340 	val p3= #1(mv_division(denom,p2,LEX_));

  2341 	val p1var=get_vars(p1');

  2342     in

  2343 	if length(xs)>0 then

  2344 	    if p3=[(1,mv_null2(var))] then

  2345 		(

  2346 		 Const ("op +",[HOLogic.realT,HOLogic.realT]--->HOLogic.realT) $

  2347 		 poly2term(the(term2poly p1' p1var),p1var) $

  2348 		 com_den2(xs,denom,den,var)

  2349 		)

  2350 	    else

  2351 		(

  2352 		 Const ("op +",[HOLogic.realT,HOLogic.realT]--->HOLogic.realT) $

  2353 		 (

  2354 		   let

  2355 		       val p3'=poly2term(p3,var);

  2356 		       val vars= (((map free2str) o vars) p1') union (((map free2str) o vars) p3');

  2357 		   in

  2358 		       poly2term(sort (mv_geq LEX_) (mv_mul(the(term2poly p1' vars) ,the(term2poly p3' vars),LEX_)),vars)

  2359 		   end

  2360 		  ) $

  2361 		 com_den2(xs,denom,den,var)

  2362 		)

  2363 	else

  2364 	    if p3=[(1,mv_null2(var))] then

  2365 		(

  2366 		 poly2term(the(term2poly p1' p1var),p1var)

  2367 		 )

  2368 	     else

  2369 		 (

  2370 		   let

  2371 		       val p3'=poly2term(p3,var);

  2372 		       val vars= (((map free2str) o vars) p1') union (((map free2str) o vars) p3');

  2373 		   in

  2374 		       poly2term(sort (mv_geq LEX_) (mv_mul(the(term2poly p1' vars) ,the(term2poly p3' vars),LEX_)),vars)

  2375 		   end

  2376 		  )

  2377     end;

  2379 (* WN0210???SK brauch ma des überhaupt *)

  2380 fun com_den_exp2(x::xs,denom,den,var)=

  2381     let

  2382 	val (t as Const ("HOL.divide",_) $ p1' $ p2')=x;

  2383 	val p2= sort (mv_geq LEX_) (the(expanded2poly p2' var));

  2384 	val p3= #1(mv_division(denom,p2,LEX_));

  2385 	val p1var=get_vars p1';

  2386     in

  2387 	if length(xs)>0 then

  2388 	    if p3=[(1,mv_null2(var))] then

  2389 		(

  2390 		 Const ("op +",[HOLogic.realT,HOLogic.realT]--->HOLogic.realT) $

  2391 		 poly2expanded(the (expanded2poly p1' p1var),p1var) $

  2392 		 com_den_exp2(xs,denom,den,var)

  2393 		)

  2394 	    else

  2395 		(

  2396 		 Const ("op +",[HOLogic.realT,HOLogic.realT]--->HOLogic.realT) $

  2397 		 (

  2398 		   let

  2399 		       val p3'=poly2expanded(p3,var);

  2400 		       val vars= (((map free2str) o vars) p1') union (((map free2str) o vars) p3');

  2401 		   in

  2402 		       poly2expanded(sort (mv_geq LEX_) (mv_mul(the(expanded2poly p1' vars) ,the(expanded2poly p3' vars),LEX_)),vars)

  2403 		   end

  2404 		  ) $

  2405 		 com_den_exp2(xs,denom,den,var)

  2406 		)

  2407 	else

  2408 	    if p3=[(1,mv_null2(var))] then

  2409 		(

  2410 		 poly2expanded(the (expanded2poly p1' p1var),p1var)

  2411 		 )

  2412 	     else

  2413 		 (

  2414 		   let

  2415 		       val p3'=poly2expanded(p3,var);

  2416 		       val vars= (((map free2str) o vars) p1') union (((map free2str) o vars) p3');

  2417 		   in

  2418 		       poly2expanded(sort (mv_geq LEX_) (mv_mul(the(expanded2poly p1' vars) ,the(expanded2poly p3' vars),LEX_)),vars)

  2419 		   end

  2420 		  )

  2421     end;

  2422 ---------------------------------------------------------*)

  2425 (*. searches for an element y of a list ys, which has an gcd not 1 with x .*)

  2426 fun exists_gcd (x,[]) = false

  2427   | exists_gcd (x,y::ys) = if mv_gcd x y = [(1,mv_null2(#2(hd(x))))] then  exists_gcd (x,ys)

  2428 			   else true;

  2430 (*. divides each element of the list xs with y .*)

  2431 fun list_div ([],y) = []

  2432   | list_div (x::xs,y) =

  2433     let

  2434 	val (d,r)=mv_division(x,y,LEX_);

  2435     in

  2436 	if r=[] then

  2437 	    d::list_div(xs,y)

  2438 	else x::list_div(xs,y)

  2439     end;

  2441 (*. checks if x is in the list ys .*)

  2442 fun in_list (x,[]) = false

  2443   | in_list (x,y::ys) = if x=y then true

  2444 			else in_list(x,ys);

  2446 (*. deletes all equal elements of the list xs .*)

  2447 fun kill_equal [] = []

  2448   | kill_equal (x::xs) = if in_list(x,xs) orelse x=[(1,mv_null2(#2(hd(x))))] then kill_equal(xs)

  2449 			 else x::kill_equal(xs);

  2451 (*. searches for new factors .*)

  2452 fun new_factors [] = []

  2453   | new_factors (list:mv_poly list):mv_poly list =

  2454     let

  2455 	val l = kill_equal list;

  2456 	val len = length(l);

  2457     in

  2458 	if len>=2 then

  2459 	    (

  2460 	     let

  2461 		 val x::y::xs=l;

  2462 		 val gcd=mv_gcd x y;

  2463 	     in

  2464 		 if gcd=[(1,mv_null2(#2(hd(x))))] then

  2465 		     (

  2466 		      if exists_gcd(x,xs) then new_factors (y::xs @ [x])

  2467 		      else x::new_factors(y::xs)

  2468 	             )

  2469 		 else gcd::new_factors(kill_equal(list_div(x::y::xs,gcd)))

  2470 	     end

  2471 	     )

  2472 	else

  2473 	    if len=1 then [hd(l)]

  2474 	    else []

  2475     end;

  2477 (*. gets the factors of a list .*)

  2478 fun get_factors x = new_factors x;

  2480 (*. multiplies the elements of the list .*)

  2481 fun multi_list [] = []

  2482   | multi_list (x::xs) = if xs=[] then x

  2483 			 else mv_mul(x,multi_list xs,LEX_);

  2485 (*. makes a term out of the elements of the list (polynomial representation) .*)

  2486 fun make_term ([],vars) = Free(str_of_int 0,HOLogic.realT)

  2487   | make_term ((x::xs),vars) = if length(xs)=0 then poly2term(sort (mv_geq LEX_) (x),vars)

  2488 			       else

  2489 				   (

  2490 				    Const ("op *",[HOLogic.realT,HOLogic.realT]--->HOLogic.realT) $

  2491 				    poly2term(sort (mv_geq LEX_) (x),vars) $

  2492 				    make_term(xs,vars)

  2493 				    );

  2495 (*. factorizes the denominator (polynomial representation) .*)

  2496 fun factorize_den (l,den,vars) =

  2497     let

  2498 	val factor_list=kill_equal( (get_factors l));

  2499 	val mlist=multi_list(factor_list);

  2500 	val (last,rest)=mv_division(den,multi_list(factor_list),LEX_);

  2501     in

  2502 	if rest=[] then

  2503 	    (

  2504 	     if last=[(1,mv_null2(vars))] then make_term(factor_list,vars)

  2505 	     else make_term(last::factor_list,vars)

  2506 	     )

  2507 	else raise error ("RATIONALS_FACTORIZE_DEN_EXCEPTION: Invalid factor by division")

  2508     end;

  2510 (*. makes a term out of the elements of the list (expanded polynomial representation) .*)

  2511 fun make_exp ([],vars) = Free(str_of_int 0,HOLogic.realT)

  2512   | make_exp ((x::xs),vars) = if length(xs)=0 then poly2expanded(sort (mv_geq LEX_) (x),vars)

  2513 			       else

  2514 				   (

  2515 				    Const ("op *",[HOLogic.realT,HOLogic.realT]--->HOLogic.realT) $

  2516 				    poly2expanded(sort (mv_geq LEX_) (x),vars) $

  2517 				    make_exp(xs,vars)

  2518 				    );

  2520 (*. factorizes the denominator (expanded polynomial representation) .*)

  2521 fun factorize_den_exp (l,den,vars) =

  2522     let

  2523 	val factor_list=kill_equal( (get_factors l));

  2524 	val mlist=multi_list(factor_list);

  2525 	val (last,rest)=mv_division(den,multi_list(factor_list),LEX_);

  2526     in

  2527 	if rest=[] then

  2528 	    (

  2529 	     if last=[(1,mv_null2(vars))] then make_exp(factor_list,vars)

  2530 	     else make_exp(last::factor_list,vars)

  2531 	     )

  2532 	else raise error ("RATIONALS_FACTORIZE_DEN_EXP_EXCEPTION: Invalid factor by division")

  2533     end;

  2535 (*. calculates the common denominator of all elements of the list and multiplies .*)

  2536 (*. the nominators and denominators with the correct factor .*)

  2537 (*. (polynomial representation) .*)

  2538 fun step_add_list_of_fractions []=(Free("0",HOLogic.realT),[]:term list)

  2539   | step_add_list_of_fractions [x]= raise error ("RATIONALS_STEP_ADD_LIST_OF_FRACTIONS_EXCEPTION: Nothing to add")

  2540   | step_add_list_of_fractions (xs) =

  2541     let

  2542         val den_list=termlist2denominators (xs); (* list of denominators *)

  2543 	val (denom,var)=calc_lcm(den_list);      (* common denominator *)

  2544 	val den=factorize_den(#1(den_list),denom,var); (* faktorisierter Nenner !!! *)

  2545     in

  2546 	com_den(xs,denom,den,var)

  2547     end;

  2549 (*. calculates the common denominator of all elements of the list and multiplies .*)

  2550 (*. the nominators and denominators with the correct factor .*)

  2551 (*. (expanded polynomial representation) .*)

  2552 fun step_add_list_of_fractions_exp []  = (Free("0",HOLogic.realT),[]:term list)

  2553   | step_add_list_of_fractions_exp [x] = raise error ("RATIONALS_STEP_ADD_LIST_OF_FRACTIONS_EXP_EXCEPTION: Nothing to add")

  2554   | step_add_list_of_fractions_exp (xs)=

  2555     let

  2556         val den_list=termlist2denominators_exp (xs); (* list of denominators *)

  2557 	val (denom,var)=calc_lcm(den_list);      (* common denominator *)

  2558 	val den=factorize_den_exp(#1(den_list),denom,var); (* faktorisierter Nenner !!! *)

  2559     in

  2560 	com_den_exp(xs,denom,den,var)

  2561     end;

  2563 (* wird aktuell nicht mehr gebraucht, bei rückänderung schon

  2564 -------------------------------------------------------------

  2565 (* WN0210???SK brauch ma des überhaupt *)

  2566 fun step_add_list_of_fractions2 []=(Free("0",HOLogic.realT),[]:term list)

  2567   | step_add_list_of_fractions2 [x]=(x,[])

  2568   | step_add_list_of_fractions2 (xs) =

  2569     let

  2570         val den_list=termlist2denominators (xs); (* list of denominators *)

  2571 	val (denom,var)=calc_lcm(den_list);      (* common denominator *)

  2572 	val den=factorize_den(#1(den_list),denom,var);  (* faktorisierter Nenner !!! *)

  2573     in

  2574 	(

  2575 	 Const ("HOL.divide",[HOLogic.realT,HOLogic.realT]--->HOLogic.realT) $

  2576 	 com_den2(xs,denom, poly2term(denom,var)(*den*),var) $

  2577 	 poly2term(denom,var)

  2578 	,

  2579 	[]

  2580 	)

  2581     end;

  2583 (* WN0210???SK brauch ma des überhaupt *)

  2584 fun step_add_list_of_fractions2_exp []=(Free("0",HOLogic.realT),[]:term list)

  2585   | step_add_list_of_fractions2_exp [x]=(x,[])

  2586   | step_add_list_of_fractions2_exp (xs) =

  2587     let

  2588         val den_list=termlist2denominators_exp (xs); (* list of denominators *)

  2589 	val (denom,var)=calc_lcm(den_list);      (* common denominator *)

  2590 	val den=factorize_den_exp(#1(den_list),denom,var);  (* faktorisierter Nenner !!! *)

  2591     in

  2592 	(

  2593 	 Const ("HOL.divide",[HOLogic.realT,HOLogic.realT]--->HOLogic.realT) $

  2594 	 com_den_exp2(xs,denom, poly2term(denom,var)(*den*),var) $

  2595 	 poly2expanded(denom,var)

  2596 	,

  2597 	[]

  2598 	)

  2599     end;

  2600 ---------------------------------------------- *)

  2603 (*. converts a term, which contains severel terms seperated by +, into a list of these terms .*)

  2604 fun term2list (t as (Const("HOL.divide",_) $ _ $ _)) = [t]

  2605   | term2list (t as (Const("Atools.pow",_) $ _ $ _)) =

  2606     [Const("HOL.divide",[HOLogic.realT,HOLogic.realT]--->HOLogic.realT) $

  2607 	  t $ Free("1",HOLogic.realT)

  2608      ]

  2609   | term2list (t as (Free(_,_))) =

  2610     [Const("HOL.divide",[HOLogic.realT,HOLogic.realT]--->HOLogic.realT) $

  2611 	  t $  Free("1",HOLogic.realT)

  2612      ]

  2613   | term2list (t as (Const("op *",_) $ _ $ _)) =

  2614     [Const("HOL.divide",[HOLogic.realT,HOLogic.realT]--->HOLogic.realT) $

  2615 	  t $ Free("1",HOLogic.realT)

  2616      ]

  2617   | term2list (Const("op +",_) $ t1 $ t2) = term2list(t1) @ term2list(t2)

  2618   | term2list (Const("op -",_) $ t1 $ t2) =

  2619     raise error ("RATIONALS_TERM2LIST_EXCEPTION: - not implemented yet")

  2620   | term2list _ = raise error ("RATIONALS_TERM2LIST_EXCEPTION: invalid term");

  2622 (*.factors out the gcd of nominator and denominator:

  2623    a/b = (a' * gcd)/(b' * gcd),  a,b,gcd  are poly[2].*)

  2624 fun factout_p_  (thy:theory) t = Some (step_cancel t,[]:term list);

  2625 fun factout_ (thy:theory) t = Some (step_cancel_expanded t,[]:term list);

  2627 (*.cancels a single fraction with normalform [2]

  2628    resulting in a canceled fraction [2], see factout_ .*)

  2629 fun cancel_p_ (thy:theory) t = (*WN.2.6.03 no rewrite -> None !*)

  2630     (let val (t',asm) = direct_cancel(*_expanded ... corrected MG.21.8.03*) t

  2631      in if t = t' then None else Some (t',asm)

  2632      end) handle _ => None;

  2633 (*.the same as above with normalform [3]

  2634   val cancel_ :

  2635       theory ->        (*10.02 unused                                    *)

  2636       term -> 	       (*fraction in normalform [3]                      *)

  2637       (term * 	       (*fraction in normalform [3]                      *)

  2638        term list)      (*casual asumptions in normalform [3]             *)

  2639 	  option       (*None: the function is not applicable            *).*)

  2640 fun cancel_ (thy:theory) t = Some (direct_cancel_expanded t) handle _ => None;

  2642 (*.transforms sums of at least 2 fractions [3] to

  2643    sums with the least common multiple as nominator.*)

  2644 fun common_nominator_p_ (thy:theory) t =

  2645 ((*writeln("### common_nominator_p_ called");*)

  2646     Some (step_add_list_of_fractions(term2list(t))) handle _ => None

  2647 );

  2648 fun common_nominator_ (thy:theory) t =

  2649     Some (step_add_list_of_fractions_exp(term2list(t))) handle _ => None;

  2651 (*.add 2 or more fractions

  2652 val add_fraction_p_ :

  2653       theory ->        (*10.02 unused                                    *)

  2654       term -> 	       (*2 or more fractions with normalform [2]         *)

  2655       (term * 	       (*one fraction with normalform [2]                *)

  2656        term list)      (*casual assumptions in normalform [2] WN0210???SK  *)

  2657 	  option       (*None: the function is not applicable            *).*)

  2658 fun add_fraction_p_ (thy:theory) t =

  2659 (writeln("### add_fraction_p_ called");

  2660     (let val ts = term2list t

  2661      in if 1 < length ts

  2662 	then Some (add_list_of_fractions ts)

  2663 	else None (*raise error ("RATIONALS_ADD_EXCEPTION: nothing to add")*)

  2664      end) handle _ => None

  2665 );

  2666 (*.same as add_fraction_p_ but with normalform [3].*)

  2667 (*Some (step_add_list_of_fractions2(term2list(t))); *)

  2668 fun add_fraction_ (thy:theory) t =

  2669     if length(term2list(t))>1

  2670     then Some (add_list_of_fractions_exp(term2list(t))) handle _ => None

  2671     else (*raise error ("RATIONALS_ADD_FRACTION_EXCEPTION: nothing to add")*)

  2672 	None;

  2673 fun add_fraction_ (thy:theory) t =

  2674     (if 1 < length (term2list t)

  2675      then Some (add_list_of_fractions_exp (term2list t))

  2676      else (*raise error ("RATIONALS_ADD_FRACTION_EXCEPTION: nothing to add")*)

  2677 	 None) handle _ => None;

  2679 (*Some (step_add_list_of_fractions2_exp(term2list(t))); *)

  2681 (*. brings the term into a normal form .*)

  2682 fun norm_rational_ (thy:theory) t =

  2683     Some (add_list_of_fractions(term2list(t))) handle _ => None;

  2684 fun norm_expanded_rat_ (thy:theory) t =

  2685     Some (add_list_of_fractions_exp(term2list(t))) handle _ => None;

  2688 (*.evaluates conditions in calculate_Rational.*)

  2689 (*make local with FIXX@ME result:term *term list*)

  2690 val calc_rat_erls = prep_rls(

  2691   Rls {id = "calc_rat_erls", preconds = [], rew_ord = ("dummy_ord",dummy_ord),

  2692 	 erls = e_rls, srls = Erls, calc = [], (*asm_thm = [], *)

  2693 	 rules =

  2694 	 [Calc ("op =",eval_equal "#equal_"),

  2695 	  Calc ("Atools.is'_const",eval_const "#is_const_"),

  2696 	  Thm ("not_true",num_str not_true),

  2697 	  Thm ("not_false",num_str not_false)

  2698 	  ],

  2699 	 scr = EmptyScr});

  2702 (*.simplifies expressions with numerals;

  2703    does NOT rearrange the term by AC-rewriting; thus terms with variables

  2704    need to have constants to be commuted together respectively.*)

  2705 val calculate_Rational = prep_rls(

  2706     merge_rls "calculate_Rational"

  2707 	(Rls {id = "divide", preconds = [], rew_ord = ("dummy_ord",dummy_ord),

  2708 	      erls = calc_rat_erls, srls = Erls, (*asm_thm = [],*)

  2709 	      calc = [],

  2710 	      rules =

  2711 	      [Calc ("HOL.divide"  ,eval_cancel "#divide_"),

  2713 	       Thm ("sym_real_minus_divide_eq",

  2714 		    num_str (real_minus_divide_eq RS sym)),

  2715 	       (*SYM - ?x / ?y = - (?x / ?y)  may come from subst*)

  2717 	       Thm ("rat_add",num_str rat_add),

  2718 	       (*"[| a is_const; b is_const; c is_const; d is_const |] ==> \

  2719 		 \"a / c + b / d = (a * d) / (c * d) + (b * c ) / (d * c)"*)

  2720 	       Thm ("rat_add1",num_str rat_add1),

  2721 	       (*"[| a is_const; b is_const; c is_const |] ==> \

  2722 		 \"a / c + b / c = (a + b) / c"*)

  2723 	       Thm ("rat_add2",num_str rat_add2),

  2724 	       (*"[| ?a is_const; ?b is_const; ?c is_const |] ==> \

  2725 		 \?a / ?c + ?b = (?a + ?b * ?c) / ?c"*)

  2726 	       Thm ("rat_add3",num_str rat_add3),

  2727 	       (*"[| a is_const; b is_const; c is_const |] ==> \

  2728 		 \"a + b / c = (a * c) / c + b / c"\

  2729 		 \.... is_const to be omitted here FIXME*)

  2731 	       Thm ("rat_mult",num_str rat_mult),

  2732 	       (*a / b * (c / d) = a * c / (b * d)*)

  2733 	       Thm ("real_times_divide1_eq",num_str real_times_divide1_eq),

  2734 	       (*?x * (?y / ?z) = ?x * ?y / ?z*)

  2735 	       Thm ("real_times_divide2_eq",num_str real_times_divide2_eq),

  2736 	       (*?y / ?z * ?x = ?y * ?x / ?z*)

  2738 	       Thm ("real_divide_divide1",num_str real_divide_divide1),

  2739 	       (*"?y ~= 0 ==> ?u / ?v / (?y / ?z) = ?u / ?v * (?z / ?y)"*)

  2740 	       Thm ("real_divide_divide2_eq",num_str real_divide_divide2_eq),

  2741 	       (*"?x / ?y / ?z = ?x / (?y * ?z)"*)

  2743 	       Thm ("rat_power", num_str rat_power),

  2744 	       (*"(?a / ?b) ^^^ ?n = ?a ^^^ ?n / ?b ^^^ ?n"*)

  2746 	       Thm ("mult_cross",num_str mult_cross),

  2747 	       (*"[| b ~= 0; d ~= 0 |] ==> (a / b = c / d) = (a * d = b * c)*)

  2748 	       Thm ("mult_cross1",num_str mult_cross1),

  2749 	       (*"   b ~= 0            ==> (a / b = c    ) = (a     = b * c)*)

  2750 	       Thm ("mult_cross2",num_str mult_cross2)

  2751 	       (*"           d ~= 0    ==> (a     = c / d) = (a * d =     c)*)

  2752 	       ], scr = EmptyScr})

  2753 	calculate_Poly);

  2756 (*("is_expanded", ("Rational.is'_expanded", eval_is_expanded ""))*)

  2757 fun eval_is_expanded (thmid:string) _

  2758 		       (t as (Const("Rational.is'_expanded", _) $ arg)) thy =

  2759     if is_expanded arg

  2760     then Some (mk_thmid thmid ""

  2761 			((string_of_cterm o cterm_of (sign_of thy)) arg) "",

  2762 	       Trueprop $ (mk_equality (t, HOLogic.true_const)))

  2763     else Some (mk_thmid thmid ""

  2764 			((string_of_cterm o cterm_of (sign_of thy)) arg) "",

  2765 	       Trueprop $ (mk_equality (t, HOLogic.false_const)))

  2766   | eval_is_expanded _ _ _ _ = None;

  2768 val rational_erls =

  2769     merge_rls "rational_erls" calculate_Rational

  2770 	      (append_rls "is_expanded" Atools_erls

  2771 			  [Calc ("Rational.is'_expanded", eval_is_expanded "")

  2772 			   ]);

  2776 (*.3 'reverse-rewrite-sets' for symbolic computation on rationals:

  2777  =================================================================

  2778  A[2] 'cancel_p': .

  2779  A[3] 'cancel': .

  2780  B[2] 'common_nominator_p': transforms summands in a term [2]

  2781          to fractions with the (least) common multiple as nominator.

  2782  B[3] 'norm_rational': normalizes arbitrary algebraic terms (without

  2783          radicals and transzendental functions) to one canceled fraction,

  2784 	 nominator and denominator in polynomial form.

  2786 In order to meet isac's requirements for interactive and stepwise calculation,

  2787 each 'reverse-rewerite-set' consists of an initialization for the interpreter

  2788 state and of 4 functions, each of which employs rewriting as much as possible.

  2789 The signature of these functions are the same in each 'reverse-rewrite-set'

  2790 respectively.*)

  2792 (* ************************************************************************* *)

  2795 local(*. cancel_p

  2796 ------------------------

  2797 cancels a single fraction consisting of two (uni- or multivariate)

  2798 polynomials WN0609???SK[2] into another such a fraction; examples:

  2800 	   a^2 + -1*b^2         a + b

  2801         -------------------- = ---------

  2802 	a^2 + -2*a*b + b^2     a + -1*b

  2804         a^2    a

  2805         --- = ---

  2806          a     1

  2808 Remark: the reverse ruleset does _NOT_ work properly with other input !.*)

  2809 (*WN020824 wir werden "uberlegen, wie wir ungeeignete inputs zur"uckweisen*)

  2811 val {rules, rew_ord=(_,ro),...} =

  2812     rep_rls (assoc_rls "make_polynomial");

  2813 (*WN060829 ... make_deriv does not terminate with 1st expl above,

  2814            see rational.sml --- investigate rulesets for cancel_p ---*)

  2815 val {rules, rew_ord=(_,ro),...} =

  2816     rep_rls (assoc_rls "rev_rew_p");

  2818 val thy = Rational.thy;

  2820 (*.init_state = fn : term -> istate

  2821 initialzies the state of the script interpreter. The state is:

  2823 type rrlsstate =      (*state for reverse rewriting*)

  2824      (term *          (*the current formula*)

  2825       term *          (*the final term*)

  2826       rule list       (*'reverse rule list' (#)*)

  2827 	    list *    (*may be serveral, eg. in norm_rational*)

  2828       (rule *         (*Thm (+ Thm generated from Calc) resulting in ...*)

  2829        (term *        (*... rewrite with ...*)

  2830 	term list))   (*... assumptions*)

  2831 	  list);      (*derivation from given term to normalform

  2832 		       in reverse order with sym_thm;

  2833                        (#) could be extracted from here by (map #1)*).*)

  2834 (* val {rules, rew_ord=(_,ro),...} =

  2835        rep_rls (assoc_rls "rev_rew_p")        (*USE ALWAYS, SEE val cancel_p*);

  2836    val (thy, eval_rls, ro) =(Rational.thy, Atools_erls, ro) (*..val cancel_p*);

  2837    val t = t;

  2838    *)

  2839 fun init_state thy eval_rls ro t =

  2840     let val Some (t',_) = factout_p_ thy t

  2841         val Some (t'',asm) = cancel_p_ thy t

  2842         val der = reverse_deriv thy eval_rls rules ro None t'

  2843         val der = der @ [(Thm ("real_mult_div_cancel2",

  2844 			       num_str real_mult_div_cancel2),

  2845 			  (t'',asm))]

  2846         val rs = (distinct_Thm o (map #1)) der

  2847 	val rs = filter_out (eq_Thms ["sym_real_add_zero_left",

  2848 				      "sym_real_mult_0",

  2849 				      "sym_real_mult_1"

  2850 				      (*..insufficient,eg.make_Polynomial*)])rs

  2851     in (t,t'',[rs(*here only _ONE_ to ease locate_rule*)],der) end;

  2853 (*.locate_rule = fn : rule list -> term -> rule

  2854 		      -> (rule * (term * term list) option) list.

  2855   checks a rule R for being a cancel-rule, and if it is,

  2856   then return the list of rules (+ the terms they are rewriting to)

  2857   which need to be applied before R should be applied.

  2858   precondition: the rule is applicable to the argument-term.

  2859 arguments:

  2860   rule list: the reverse rule list

  2861   -> term  : ... to which the rule shall be applied

  2862   -> rule  : ... to be applied to term

  2863 value:

  2864   -> (rule           : a rule rewriting to ...

  2865       * (term        : ... the resulting term ...

  2866          * term list): ... with the assumptions ( //#0).

  2867       ) list         : there may be several such rules;

  2868 		       the list is empty, if the rule has nothing to do

  2869 		       with cancelation.*)

  2870 (* val () = ();

  2871    *)

  2872 fun locate_rule thy eval_rls ro [rs] t r =

  2873     if (id_of_thm r) mem (map (id_of_thm)) rs

  2874     then let val ropt =

  2875 		 rewrite_ thy ro eval_rls true (thm_of_thm r) t;

  2876 	 in case ropt of

  2877 		Some ta => [(r, ta)]

  2878 	      | None => (writeln("### locate_rule:  rewrite "^

  2879 				 (id_of_thm r)^" "^(term2str t)^" = None");

  2880 			 []) end

  2881     else (writeln("### locate_rule:  "^(id_of_thm r)^" not mem rrls");[])

  2882   | locate_rule _ _ _ _ _ _ =

  2883     raise error ("locate_rule: doesnt match rev-sets in istate");

  2885 (*.next_rule = fn : rule list -> term -> rule option

  2886   for a given term return the next rules to be done for cancelling.

  2887 arguments:

  2888   rule list     : the reverse rule list

  2889   term          : the term for which ...

  2890 value:

  2891   -> rule option: ... this rule is appropriate for cancellation;

  2892 		  there may be no such rule (if the term is canceled already.*)

  2893 (* val thy = Rational.thy;

  2894    val Rrls {rew_ord=(_,ro),...} = cancel;

  2895    val ([rs],t) = (rss,f);

  2896    next_rule thy eval_rls ro [rs] t;(*eval fun next_rule ... before!*)

  2898    val (thy, [rs]) = (Rational.thy, revsets);

  2899    val Rrls {rew_ord=(_,ro),...} = cancel;

  2900    nex [rs] t;

  2901    *)

  2902 fun next_rule thy eval_rls ro [rs] t =

  2903     let val der = make_deriv thy eval_rls rs ro None t;

  2904     in case der of

  2905 (* val (_,r,_)::_ = der;

  2906    *)

  2907 	   (_,r,_)::_ => Some r

  2908 	 | _ => None

  2909     end

  2910   | next_rule _ _ _ _ _ =

  2911     raise error ("next_rule: doesnt match rev-sets in istate");

  2913 (*.val attach_form = f : rule list -> term -> term

  2914 			 -> (rule * (term * term list)) list

  2915   checks an input term TI, if it may belong to a current cancellation, by

  2916   trying to derive it from the given term TG.

  2917 arguments:

  2918   term   : TG, the last one in the cancellation agreed upon by user + math-eng

  2919   -> term: TI, the next one input by the user

  2920 value:

  2921   -> (rule           : the rule to be applied in order to reach TI

  2922       * (term        : ... obtained by applying the rule ...

  2923          * term list): ... and the respective assumptions.

  2924       ) list         : there may be several such rules;

  2925                        the list is empty, if the users term does not belong

  2926 		       to a cancellation of the term last agreed upon.*)

  2927 fun attach_form (_:rule list list) (_:term) (_:term) = (*still missing*)

  2928     []:(rule * (term * term list)) list;

  2930 in

  2932 val cancel_p =

  2933     Rrls {id = "cancel_p", prepat=[],

  2934 	  rew_ord=("ord_make_polynomial",

  2935 		   ord_make_polynomial false Rational.thy),

  2936 	  erls = rational_erls,

  2937 	  calc = [("plus"    ,("op +"        ,eval_binop "#add_")),

  2938 		  ("times"   ,("op *"        ,eval_binop "#mult_")),

  2939 		  ("divide_" ,("HOL.divide"  ,eval_cancel "#divide_")),

  2940 		  ("power_"  ,("Atools.pow"  ,eval_binop "#power_"))],

  2941 	  (*asm_thm=[("real_mult_div_cancel2","")],*)

  2942 	  scr=Rfuns {init_state  = init_state thy Atools_erls ro,

  2943 		     normal_form = cancel_p_ thy,

  2944 		     locate_rule = locate_rule thy Atools_erls ro,

  2945 		     next_rule   = next_rule thy Atools_erls ro,

  2946 		     attach_form = attach_form}}

  2947 end;(*local*)

  2950 local(*.ad (1) 'cancel'

  2951 ------------------------------

  2952 cancels a single fraction consisting of two (uni- or multivariate)

  2953 polynomials WN0609???SK[3] into another such a fraction; examples:

  2955 	   a^2 - b^2           a + b

  2956         -------------------- = ---------

  2957 	a^2 - 2*a*b + b^2      a - *b

  2959 Remark: the reverse ruleset does _NOT_ work properly with other input !.*)

  2960 (*WN 24.8.02: wir werden "uberlegen, wie wir ungeeignete inputs zur"uckweisen*)

  2962 (*

  2963 val Some (Rls {rules=rules,rew_ord=(_,ro),...}) =

  2964     assoc'(!ruleset',"expand_binoms");

  2965 *)

  2966 val {rules=rules,rew_ord=(_,ro),...} =

  2967     rep_rls (assoc_rls "expand_binoms");

  2968 val thy = Rational.thy;

  2970 fun init_state thy eval_rls ro t =

  2971     let val Some (t',_) = factout_ thy t;

  2972         val Some (t'',asm) = cancel_ thy t;

  2973         val der = reverse_deriv thy eval_rls rules ro None t';

  2974         val der = der @ [(Thm ("real_mult_div_cancel2",

  2975 			       num_str real_mult_div_cancel2),

  2976 			  (t'',asm))]

  2977         val rs = map #1 der;

  2978     in (t,t'',[rs],der) end;

  2980 fun locate_rule thy eval_rls ro [rs] t r =

  2981     if (id_of_thm r) mem (map (id_of_thm)) rs

  2982     then let val ropt =

  2983 		 rewrite_ thy ro eval_rls true (thm_of_thm r) t;

  2984 	 in case ropt of

  2985 		Some ta => [(r, ta)]

  2986 	      | None => (writeln("### locate_rule:  rewrite "^

  2987 				 (id_of_thm r)^" "^(term2str t)^" = None");

  2988 			 []) end

  2989     else (writeln("### locate_rule:  "^(id_of_thm r)^" not mem rrls");[])

  2990   | locate_rule _ _ _ _ _ _ =

  2991     raise error ("locate_rule: doesnt match rev-sets in istate");

  2993 fun next_rule thy eval_rls ro [rs] t =

  2994     let val der = make_deriv thy eval_rls rs ro None t;

  2995     in case der of

  2996 (* val (_,r,_)::_ = der;

  2997    *)

  2998 	   (_,r,_)::_ => Some r

  2999 	 | _ => None

  3000     end

  3001   | next_rule _ _ _ _ _ =

  3002     raise error ("next_rule: doesnt match rev-sets in istate");

  3004 fun attach_form (_:rule list list) (_:term) (_:term) = (*still missing*)

  3005     []:(rule * (term * term list)) list;

  3007 val pat = (term_of o the o (parse thy)) "?r/?s";

  3008 val pre1 = (term_of o the o (parse thy)) "r is_expanded";

  3009 val pre2 = (term_of o the o (parse thy)) "s is_expanded";

  3010 val prepat = [([pre1, pre2], pat)];

  3012 in

  3015 val cancel =

  3016     Rrls {id = "cancel", prepat=prepat,

  3017 	  rew_ord=("ord_make_polynomial",

  3018 		   ord_make_polynomial false Rational.thy),

  3019 	  erls = rational_erls,

  3020 	  calc = [("plus"    ,("op +"        ,eval_binop "#add_")),

  3021 		  ("times"   ,("op *"        ,eval_binop "#mult_")),

  3022 		  ("divide_" ,("HOL.divide"  ,eval_cancel "#divide_")),

  3023 		  ("power_"  ,("Atools.pow"  ,eval_binop "#power_"))],

  3024 	  scr=Rfuns {init_state  = init_state thy Atools_erls ro,

  3025 		     normal_form = cancel_ thy,

  3026 		     locate_rule = locate_rule thy Atools_erls ro,

  3027 		     next_rule   = next_rule thy Atools_erls ro,

  3028 		     attach_form = attach_form}}

  3029 end;(*local*)

  3033 local(*.ad [2] 'common_nominator_p'

  3034 ---------------------------------

  3035 FIXME Beschreibung .*)

  3038 val {rules=rules,rew_ord=(_,ro),...} =

  3039     rep_rls (assoc_rls "make_polynomial");

  3040 (*WN060829 ... make_deriv does not terminate with 1st expl above,

  3041            see rational.sml --- investigate rulesets for cancel_p ---*)

  3042 val {rules, rew_ord=(_,ro),...} =

  3043     rep_rls (assoc_rls "rev_rew_p");

  3044 val thy = Rational.thy;

  3047 (*.common_nominator_p_ = fn : theory -> term -> (term * term list) option

  3048   as defined above*)

  3050 (*.init_state = fn : term -> istate

  3051 initialzies the state of the interactive interpreter. The state is:

  3053 type rrlsstate =      (*state for reverse rewriting*)

  3054      (term *          (*the current formula*)

  3055       term *          (*the final term*)

  3056       rule list       (*'reverse rule list' (#)*)

  3057 	    list *    (*may be serveral, eg. in norm_rational*)

  3058       (rule *         (*Thm (+ Thm generated from Calc) resulting in ...*)

  3059        (term *        (*... rewrite with ...*)

  3060 	term list))   (*... assumptions*)

  3061 	  list);      (*derivation from given term to normalform

  3062 		       in reverse order with sym_thm;

  3063                        (#) could be extracted from here by (map #1)*).*)

  3064 fun init_state thy eval_rls ro t =

  3065     let val Some (t',_) = common_nominator_p_ thy t;

  3066         val Some (t'',asm) = add_fraction_p_ thy t;

  3067         val der = reverse_deriv thy eval_rls rules ro None t';

  3068         val der = der @ [(Thm ("real_mult_div_cancel2",

  3069 			       num_str real_mult_div_cancel2),

  3070 			  (t'',asm))]

  3071         val rs = (distinct_Thm o (map #1)) der;

  3072 	val rs = filter_out (eq_Thms ["sym_real_add_zero_left",

  3073 				      "sym_real_mult_0",

  3074 				      "sym_real_mult_1"]) rs;

  3075     in (t,t'',[rs(*here only _ONE_*)],der) end;

  3077 (* use"knowledge/Rational.ML";

  3078    *)

  3080 (*.locate_rule = fn : rule list -> term -> rule

  3081 		      -> (rule * (term * term list) option) list.

  3082   checks a rule R for being a cancel-rule, and if it is,

  3083   then return the list of rules (+ the terms they are rewriting to)

  3084   which need to be applied before R should be applied.

  3085   precondition: the rule is applicable to the argument-term.

  3086 arguments:

  3087   rule list: the reverse rule list

  3088   -> term  : ... to which the rule shall be applied

  3089   -> rule  : ... to be applied to term

  3090 value:

  3091   -> (rule           : a rule rewriting to ...

  3092       * (term        : ... the resulting term ...

  3093          * term list): ... with the assumptions ( //#0).

  3094       ) list         : there may be several such rules;

  3095 		       the list is empty, if the rule has nothing to do

  3096 		       with cancelation.*)

  3097 (* val () = ();

  3098    *)

  3099 fun locate_rule thy eval_rls ro [rs] t r =

  3100     if (id_of_thm r) mem (map (id_of_thm)) rs

  3101     then let val ropt =

  3102 		 rewrite_ thy ro eval_rls true (thm_of_thm r) t;

  3103 	 in case ropt of

  3104 		Some ta => [(r, ta)]

  3105 	      | None => (writeln("### locate_rule:  rewrite "^

  3106 				 (id_of_thm r)^" "^(term2str t)^" = None");

  3107 			 []) end

  3108     else (writeln("### locate_rule:  "^(id_of_thm r)^" not mem rrls");[])

  3109   | locate_rule _ _ _ _ _ _ =

  3110     raise error ("locate_rule: doesnt match rev-sets in istate");

  3112 (*.next_rule = fn : rule list -> term -> rule option

  3113   for a given term return the next rules to be done for cancelling.

  3114 arguments:

  3115   rule list     : the reverse rule list

  3116   term          : the term for which ...

  3117 value:

  3118   -> rule option: ... this rule is appropriate for cancellation;

  3119 		  there may be no such rule (if the term is canceled already.*)

  3120 (* val thy = Rational.thy;

  3121    val Rrls {rew_ord=(_,ro),...} = cancel;

  3122    val ([rs],t) = (rss,f);

  3123    next_rule thy eval_rls ro [rs] t;(*eval fun next_rule ... before!*)

  3125    val (thy, [rs]) = (Rational.thy, revsets);

  3126    val Rrls {rew_ord=(_,ro),...} = cancel;

  3127    nex [rs] t;

  3128    *)

  3129 fun next_rule thy eval_rls ro [rs] t =

  3130     let val der = make_deriv thy eval_rls rs ro None t;

  3131     in case der of

  3132 (* val (_,r,_)::_ = der;

  3133    *)

  3134 	   (_,r,_)::_ => Some r

  3135 	 | _ => None

  3136     end

  3137   | next_rule _ _ _ _ _ =

  3138     raise error ("next_rule: doesnt match rev-sets in istate");

  3140 (*.val attach_form = f : rule list -> term -> term

  3141 			 -> (rule * (term * term list)) list

  3142   checks an input term TI, if it may belong to a current cancellation, by

  3143   trying to derive it from the given term TG.

  3144 arguments:

  3145   term   : TG, the last one in the cancellation agreed upon by user + math-eng

  3146   -> term: TI, the next one input by the user

  3147 value:

  3148   -> (rule           : the rule to be applied in order to reach TI

  3149       * (term        : ... obtained by applying the rule ...

  3150          * term list): ... and the respective assumptions.

  3151       ) list         : there may be several such rules;

  3152                        the list is empty, if the users term does not belong

  3153 		       to a cancellation of the term last agreed upon.*)

  3154 fun attach_form (_:rule list list) (_:term) (_:term) = (*still missing*)

  3155     []:(rule * (term * term list)) list;

  3157 val pat0 = (term_of o the o (parse thy)) "?r/?s+?u/?v";

  3158 val pat1 = (term_of o the o (parse thy)) "?r/?s+?u   ";

  3159 val pat2 = (term_of o the o (parse thy)) "?r   +?u/?v";

  3160 val prepat = [([HOLogic.true_const], pat0),

  3161 	      ([HOLogic.true_const], pat1),

  3162 	      ([HOLogic.true_const], pat2)];

  3164 in

  3166 (*11.02 schnelle L"osung f"ur RL: Bruch auch gek"urzt;

  3167   besser w"are: auf 1 gemeinsamen Bruchstrich, Nenner und Z"ahler unvereinfacht

  3168   dh. wie common_nominator_p_, aber auf 1 Bruchstrich*)

  3169 val common_nominator_p =

  3170     Rrls {id = "common_nominator_p", prepat=prepat,

  3171 	  rew_ord=("ord_make_polynomial",

  3172 		   ord_make_polynomial false Rational.thy),

  3173 	  erls = rational_erls,

  3174 	  calc = [("plus"    ,("op +"        ,eval_binop "#add_")),

  3175 		  ("times"   ,("op *"        ,eval_binop "#mult_")),

  3176 		  ("divide_" ,("HOL.divide"  ,eval_cancel "#divide_")),

  3177 		  ("power_"  ,("Atools.pow"  ,eval_binop "#power_"))],

  3178 	  scr=Rfuns {init_state  = init_state thy Atools_erls ro,

  3179 		     normal_form = add_fraction_p_ thy,(*FIXME.WN0211*)

  3180 		     locate_rule = locate_rule thy Atools_erls ro,

  3181 		     next_rule   = next_rule thy Atools_erls ro,

  3182 		     attach_form = attach_form}}

  3183 end;(*local*)

  3186 local(*.ad [2] 'common_nominator'

  3187 ---------------------------------

  3188 FIXME Beschreibung .*)

  3191 val {rules=rules,rew_ord=(_,ro),...} =

  3192     rep_rls (assoc_rls "make_polynomial");

  3193 val thy = Rational.thy;

  3196 (*.common_nominator_ = fn : theory -> term -> (term * term list) option

  3197   as defined above*)

  3199 (*.init_state = fn : term -> istate

  3200 initialzies the state of the interactive interpreter. The state is:

  3202 type rrlsstate =      (*state for reverse rewriting*)

  3203      (term *          (*the current formula*)

  3204       term *          (*the final term*)

  3205       rule list       (*'reverse rule list' (#)*)

  3206 	    list *    (*may be serveral, eg. in norm_rational*)

  3207       (rule *         (*Thm (+ Thm generated from Calc) resulting in ...*)

  3208        (term *        (*... rewrite with ...*)

  3209 	term list))   (*... assumptions*)

  3210 	  list);      (*derivation from given term to normalform

  3211 		       in reverse order with sym_thm;

  3212                        (#) could be extracted from here by (map #1)*).*)

  3213 fun init_state thy eval_rls ro t =

  3214     let val Some (t',_) = common_nominator_ thy t;

  3215         val Some (t'',asm) = add_fraction_ thy t;

  3216         val der = reverse_deriv thy eval_rls rules ro None t';

  3217         val der = der @ [(Thm ("real_mult_div_cancel2",

  3218 			       num_str real_mult_div_cancel2),

  3219 			  (t'',asm))]

  3220         val rs = (distinct_Thm o (map #1)) der;

  3221 	val rs = filter_out (eq_Thms ["sym_real_add_zero_left",

  3222 				      "sym_real_mult_0",

  3223 				      "sym_real_mult_1"]) rs;

  3224     in (t,t'',[rs(*here only _ONE_*)],der) end;

  3226 (* use"knowledge/Rational.ML";

  3227    *)

  3229 (*.locate_rule = fn : rule list -> term -> rule

  3230 		      -> (rule * (term * term list) option) list.

  3231   checks a rule R for being a cancel-rule, and if it is,

  3232   then return the list of rules (+ the terms they are rewriting to)

  3233   which need to be applied before R should be applied.

  3234   precondition: the rule is applicable to the argument-term.

  3235 arguments:

  3236   rule list: the reverse rule list

  3237   -> term  : ... to which the rule shall be applied

  3238   -> rule  : ... to be applied to term

  3239 value:

  3240   -> (rule           : a rule rewriting to ...

  3241       * (term        : ... the resulting term ...

  3242          * term list): ... with the assumptions ( //#0).

  3243       ) list         : there may be several such rules;

  3244 		       the list is empty, if the rule has nothing to do

  3245 		       with cancelation.*)

  3246 (* val () = ();

  3247    *)

  3248 fun locate_rule thy eval_rls ro [rs] t r =

  3249     if (id_of_thm r) mem (map (id_of_thm)) rs

  3250     then let val ropt =

  3251 		 rewrite_ thy ro eval_rls true (thm_of_thm r) t;

  3252 	 in case ropt of

  3253 		Some ta => [(r, ta)]

  3254 	      | None => (writeln("### locate_rule:  rewrite "^

  3255 				 (id_of_thm r)^" "^(term2str t)^" = None");

  3256 			 []) end

  3257     else (writeln("### locate_rule:  "^(id_of_thm r)^" not mem rrls");[])

  3258   | locate_rule _ _ _ _ _ _ =

  3259     raise error ("locate_rule: doesnt match rev-sets in istate");

  3261 (*.next_rule = fn : rule list -> term -> rule option

  3262   for a given term return the next rules to be done for cancelling.

  3263 arguments:

  3264   rule list     : the reverse rule list

  3265   term          : the term for which ...

  3266 value:

  3267   -> rule option: ... this rule is appropriate for cancellation;

  3268 		  there may be no such rule (if the term is canceled already.*)

  3269 (* val thy = Rational.thy;

  3270    val Rrls {rew_ord=(_,ro),...} = cancel;

  3271    val ([rs],t) = (rss,f);

  3272    next_rule thy eval_rls ro [rs] t;(*eval fun next_rule ... before!*)

  3274    val (thy, [rs]) = (Rational.thy, revsets);

  3275    val Rrls {rew_ord=(_,ro),...} = cancel_p;

  3276    nex [rs] t;

  3277    *)

  3278 fun next_rule thy eval_rls ro [rs] t =

  3279     let val der = make_deriv thy eval_rls rs ro None t;

  3280     in case der of

  3281 (* val (_,r,_)::_ = der;

  3282    *)

  3283 	   (_,r,_)::_ => Some r

  3284 	 | _ => None

  3285     end

  3286   | next_rule _ _ _ _ _ =

  3287     raise error ("next_rule: doesnt match rev-sets in istate");

  3289 (*.val attach_form = f : rule list -> term -> term

  3290 			 -> (rule * (term * term list)) list

  3291   checks an input term TI, if it may belong to a current cancellation, by

  3292   trying to derive it from the given term TG.

  3293 arguments:

  3294   term   : TG, the last one in the cancellation agreed upon by user + math-eng

  3295   -> term: TI, the next one input by the user

  3296 value:

  3297   -> (rule           : the rule to be applied in order to reach TI

  3298       * (term        : ... obtained by applying the rule ...

  3299          * term list): ... and the respective assumptions.

  3300       ) list         : there may be several such rules;

  3301                        the list is empty, if the users term does not belong

  3302 		       to a cancellation of the term last agreed upon.*)

  3303 fun attach_form (_:rule list list) (_:term) (_:term) = (*still missing*)

  3304     []:(rule * (term * term list)) list;

  3306 val pat0 =  (term_of o the o (parse thy)) "?r/?s+?u/?v";

  3307 val pat01 = (term_of o the o (parse thy)) "?r/?s-?u/?v";

  3308 val pat1 =  (term_of o the o (parse thy)) "?r/?s+?u   ";

  3309 val pat11 = (term_of o the o (parse thy)) "?r/?s-?u   ";

  3310 val pat2 =  (term_of o the o (parse thy)) "?r   +?u/?v";

  3311 val pat21 = (term_of o the o (parse thy)) "?r   -?u/?v";

  3312 val prepat = [([HOLogic.true_const], pat0),

  3313 	      ([HOLogic.true_const], pat01),

  3314 	      ([HOLogic.true_const], pat1),

  3315 	      ([HOLogic.true_const], pat11),

  3316 	      ([HOLogic.true_const], pat2),

  3317 	      ([HOLogic.true_const], pat21)];

  3320 in

  3322 val common_nominator =

  3323     Rrls {id = "common_nominator", prepat=prepat,

  3324 	  rew_ord=("ord_make_polynomial",

  3325 		   ord_make_polynomial false Rational.thy),

  3326 	  erls = rational_erls,

  3327 	  calc = [("plus"    ,("op +"        ,eval_binop "#add_")),

  3328 		  ("times"   ,("op *"        ,eval_binop "#mult_")),

  3329 		  ("divide_" ,("HOL.divide"  ,eval_cancel "#divide_")),

  3330 		  ("power_"  ,("Atools.pow"  ,eval_binop "#power_"))],

  3331 	  (*asm_thm=[("real_mult_div_cancel2","")],*)

  3332 	  scr=Rfuns {init_state  = init_state thy Atools_erls ro,

  3333 		     normal_form = add_fraction_ (*NOT common_nominator_*) thy,

  3334 		     locate_rule = locate_rule thy Atools_erls ro,

  3335 		     next_rule   = next_rule thy Atools_erls ro,

  3336 		     attach_form = attach_form}}

  3338 end;(*local*)

  3341 (*##*)

  3342 end;(*struct*)

  3344 open RationalI;

  3345 (*##*)

  3347 (*.the expression contains + - * ^ / only ?.*)

  3348 fun is_ratpolyexp (Free _) = true

  3349   | is_ratpolyexp (Const ("op +",_) $ Free _ $ Free _) = true

  3350   | is_ratpolyexp (Const ("op -",_) $ Free _ $ Free _) = true

  3351   | is_ratpolyexp (Const ("op *",_) $ Free _ $ Free _) = true

  3352   | is_ratpolyexp (Const ("Atools.pow",_) $ Free _ $ Free _) = true

  3353   | is_ratpolyexp (Const ("HOL.divide",_) $ Free _ $ Free _) = true

  3354   | is_ratpolyexp (Const ("op +",_) $ t1 $ t2) =

  3355                ((is_ratpolyexp t1) andalso (is_ratpolyexp t2))

  3356   | is_ratpolyexp (Const ("op -",_) $ t1 $ t2) =

  3357                ((is_ratpolyexp t1) andalso (is_ratpolyexp t2))

  3358   | is_ratpolyexp (Const ("op *",_) $ t1 $ t2) =

  3359                ((is_ratpolyexp t1) andalso (is_ratpolyexp t2))

  3360   | is_ratpolyexp (Const ("Atools.pow",_) $ t1 $ t2) =

  3361                ((is_ratpolyexp t1) andalso (is_ratpolyexp t2))

  3362   | is_ratpolyexp (Const ("HOL.divide",_) $ t1 $ t2) =

  3363                ((is_ratpolyexp t1) andalso (is_ratpolyexp t2))

  3364   | is_ratpolyexp _ = false;

  3366 (*("is_ratpolyexp", ("Rational.is'_ratpolyexp", eval_is_ratpolyexp ""))*)

  3367 fun eval_is_ratpolyexp (thmid:string) _

  3368 		       (t as (Const("Rational.is'_ratpolyexp", _) $ arg)) thy =

  3369     if is_ratpolyexp arg

  3370     then Some (mk_thmid thmid ""

  3371 			((string_of_cterm o cterm_of (sign_of thy)) arg) "",

  3372 	       Trueprop $ (mk_equality (t, HOLogic.true_const)))

  3373     else Some (mk_thmid thmid ""

  3374 			((string_of_cterm o cterm_of (sign_of thy)) arg) "",

  3375 	       Trueprop $ (mk_equality (t, HOLogic.false_const)))

  3376   | eval_is_ratpolyexp _ _ _ _ = None;

  3380 (*-------------------18.3.03 --> struct <-----------vvv--*)

  3381 val add_fractions_p = common_nominator_p; (*FIXXXME:eilig f"ur norm_Rational*)

  3383 (*.discard binary minus, shift unary minus into -1*;

  3384    unary minus before numerals are put into the numeral by parsing;

  3385    contains absolute minimum of thms for context in norm_Rational .*)

  3386 val discard_minus = prep_rls(

  3387   Rls {id = "discard_minus", preconds = [], rew_ord = ("dummy_ord",dummy_ord),

  3388       erls = e_rls, srls = Erls, calc = [], (*asm_thm = [],*)

  3389       rules = [Thm ("real_diff_minus", num_str real_diff_minus),

  3390 	       (*"a - b = a + -1 * b"*)

  3391 	       Thm ("sym_real_mult_minus1",num_str (real_mult_minus1 RS sym))

  3392 	       (*- ?z = "-1 * ?z"*)

  3393 	       ],

  3394       scr = Script ((term_of o the o (parse thy))

  3395       "empty_script")

  3396       }):rls;

  3397 (*erls for calculate_Rational; make local with FIXX@ME result:term *term list*)

  3398 val powers_erls = prep_rls(

  3399   Rls {id = "powers_erls", preconds = [], rew_ord = ("dummy_ord",dummy_ord),

  3400       erls = e_rls, srls = Erls, calc = [], (*asm_thm = [],*)

  3401       rules = [Calc ("Atools.is'_atom",eval_is_atom "#is_atom_"),

  3402 	       Calc ("Atools.is'_even",eval_is_even "#is_even_"),

  3403 	       Calc ("op <",eval_equ "#less_"),

  3404 	       Thm ("not_false", not_false),

  3405 	       Thm ("not_true", not_true),

  3406 	       Calc ("op +",eval_binop "#add_")

  3407 	       ],

  3408       scr = Script ((term_of o the o (parse thy))

  3409       "empty_script")

  3410       }:rls);

  3411 (*.all powers over + distributed; atoms over * collected, other distributed

  3412    contains absolute minimum of thms for context in norm_Rational .*)

  3413 val powers = prep_rls(

  3414   Rls {id = "powers", preconds = [], rew_ord = ("dummy_ord",dummy_ord),

  3415       erls = powers_erls, srls = Erls, calc = [], (*asm_thm = [],*)

  3416       rules = [Thm ("realpow_multI", num_str realpow_multI),

  3417 	       (*"(r * s) ^^^ n = r ^^^ n * s ^^^ n"*)

  3418 	       Thm ("realpow_pow",num_str realpow_pow),

  3419 	       (*"(a ^^^ b) ^^^ c = a ^^^ (b * c)"*)

  3420 	       Thm ("realpow_oneI",num_str realpow_oneI),

  3421 	       (*"r ^^^ 1 = r"*)

  3422 	       Thm ("realpow_minus_even",num_str realpow_minus_even),

  3423 	       (*"n is_even ==> (- r) ^^^ n = r ^^^ n" ?-->discard_minus?*)

  3424 	       Thm ("realpow_minus_odd",num_str realpow_minus_odd),

  3425 	       (*"Not (n is_even) ==> (- r) ^^^ n = -1 * r ^^^ n"*)

  3427 	       (*----- collect atoms over * -----*)

  3428 	       Thm ("realpow_two_atom",num_str realpow_two_atom),

  3429 	       (*"r is_atom ==> r * r = r ^^^ 2"*)

  3430 	       Thm ("realpow_plus_1",num_str realpow_plus_1),

  3431 	       (*"r is_atom ==> r * r ^^^ n = r ^^^ (n + 1)"*)

  3432 	       Thm ("realpow_addI_atom",num_str realpow_addI_atom),

  3433 	       (*"r is_atom ==> r ^^^ n * r ^^^ m = r ^^^ (n + m)"*)

  3435 	       (*----- distribute none-atoms -----*)

  3436 	       Thm ("realpow_def_atom",num_str realpow_def_atom),

  3437 	       (*"[| 1 < n; not(r is_atom) |]==>r ^^^ n = r * r ^^^ (n + -1)"*)

  3438 	       Thm ("realpow_eq_oneI",num_str realpow_eq_oneI),

  3439 	       (*"1 ^^^ n = 1"*)

  3440 	       Calc ("op +",eval_binop "#add_")

  3441 	       ],

  3442       scr = Script ((term_of o the o (parse thy))

  3443       "empty_script")

  3444       }:rls);

  3445 (*.contains absolute minimum of thms for context in norm_Rational.*)

  3446 val rat_mult_divide = prep_rls(

  3447   Rls {id = "rat_mult_divide", preconds = [],

  3448        rew_ord = ("dummy_ord",dummy_ord),

  3449       erls = e_rls, srls = Erls, calc = [], (*asm_thm = [],*)

  3450       rules = [Thm ("rat_mult",num_str rat_mult),

  3451 	       (*(1)"?a / ?b * (?c / ?d) = ?a * ?c / (?b * ?d)"*)

  3452 	       Thm ("real_times_divide1_eq",num_str real_times_divide1_eq),

  3453 	       (*(2)"?a * (?c / ?d) = ?a * ?c / ?d" must be [2],

  3454 	       otherwise inv.to a / b / c = ...*)

  3455 	       Thm ("real_times_divide2_eq",num_str real_times_divide2_eq),

  3456 	       (*"?a / ?b * ?c = ?a * ?c / ?b" order weights x^^^n too much

  3457 		     and does not commute a / b * c ^^^ 2 !*)

  3459 	       Thm ("real_divide_divide1_eq", real_divide_divide1_eq),

  3460 	       (*"?x / (?y / ?z) = ?x * ?z / ?y"*)

  3461 	       Thm ("real_divide_divide2_eq", real_divide_divide2_eq),

  3462 	       (*"?x / ?y / ?z = ?x / (?y * ?z)"*)

  3463 	       Calc ("HOL.divide"  ,eval_cancel "#divide_")

  3464 	       ],

  3465       scr = Script ((term_of o the o (parse thy)) "empty_script")

  3466       }:rls);

  3467 (*.contains absolute minimum of thms for context in norm_Rational.*)

  3468 val reduce_0_1_2 = prep_rls(

  3469   Rls{id = "reduce_0_1_2", preconds = [], rew_ord = ("dummy_ord", dummy_ord),

  3470       erls = e_rls,srls = Erls,calc = [],(*asm_thm = [],*)

  3471       rules = [(*Thm ("real_divide_1",num_str real_divide_1),

  3472 		 "?x / 1 = ?x" unnecess.for normalform*)

  3473 	       Thm ("real_mult_1",num_str real_mult_1),

  3474 	       (*"1 * z = z"*)

  3475 	       (*Thm ("real_mult_minus1",num_str real_mult_minus1),

  3476 	       "-1 * z = - z"*)

  3477 	       (*Thm ("real_minus_mult_cancel",num_str real_minus_mult_cancel),

  3478 	       "- ?x * - ?y = ?x * ?y"*)

  3480 	       Thm ("real_mult_0",num_str real_mult_0),

  3481 	       (*"0 * z = 0"*)

  3482 	       Thm ("real_add_zero_left",num_str real_add_zero_left),

  3483 	       (*"0 + z = z"*)

  3484 	       (*Thm ("real_add_minus",num_str real_add_minus),

  3485 	       "?z + - ?z = 0"*)

  3487 	       Thm ("sym_real_mult_2",num_str (real_mult_2 RS sym)),

  3488 	       (*"z1 + z1 = 2 * z1"*)

  3489 	       Thm ("real_mult_2_assoc",num_str real_mult_2_assoc),

  3490 	       (*"z1 + (z1 + k) = 2 * z1 + k"*)

  3492 	       Thm ("real_0_divide",num_str real_0_divide)

  3493 	       (*"0 / ?x = 0"*)

  3494 	       ], scr = EmptyScr}:rls);

  3496 (*erls for calculate_Rational;

  3497   make local with FIXX@ME result:term *term list WN0609???SKMG*)

  3498 val norm_rat_erls = prep_rls(

  3499   Rls {id = "norm_rat_erls", preconds = [], rew_ord = ("dummy_ord",dummy_ord),

  3500       erls = e_rls, srls = Erls, calc = [], (*asm_thm = [],*)

  3501       rules = [Calc ("Atools.is'_const",eval_const "#is_const_")

  3502 	       ],

  3503       scr = Script ((term_of o the o (parse thy))

  3504       "empty_script")

  3505       }:rls);

  3506 (*.consists of rls containing the absolute minimum of thms.*)

  3507 (*040209: this version has been used by RL for his equations,

  3508 which is now replaced by MGs version below

  3509 vvv OLD VERSION !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!*)

  3510 val norm_Rational = prep_rls(

  3511   Rls {id = "norm_Rational", preconds = [], rew_ord = ("dummy_ord",dummy_ord),

  3512       erls = norm_rat_erls, srls = Erls, calc = [], (*asm_thm = [],*)

  3513       rules = [(*sequence given by operator precedence*)

  3514 	       Rls_ discard_minus,

  3515 	       Rls_ powers,

  3516 	       Rls_ rat_mult_divide,

  3517 	       Rls_ expand,

  3518 	       Rls_ reduce_0_1_2,

  3519 	       (*^^^^^^^^^ from RL -- not the latest one vvvvvvvvv*)

  3520 	       Rls_ order_add_mult,

  3521 	       Rls_ collect_numerals,

  3522 	       Rls_ add_fractions_p,

  3523 	       Rls_ cancel_p

  3524 	       ],

  3525       scr = Script ((term_of o the o (parse thy))

  3526       "empty_script")

  3527       }:rls);

  3528 val norm_Rational_parenthesized = prep_rls(

  3529   Seq {id = "norm_Rational_parenthesized", preconds = []:term list,

  3530        rew_ord = ("dummy_ord", dummy_ord),

  3531       erls = Atools_erls, srls = Erls,

  3532       calc = [], (*asm_thm = [],*)

  3533       rules = [Rls_  norm_Rational, (*from RL -- not the latest one*)

  3534 	       Rls_ discard_parentheses

  3535 	       ],

  3536       scr = EmptyScr

  3537       }:rls);

  3540 (*-------------------18.3.03 --> struct <-----------^^^--*)

  3544 theory' := overwritel (!theory', [("Rational.thy",Rational.thy)]);

  3547 (*WN030318???SK: simplifies all but cancel and common_nominator*)

  3548 val simplify_rational =

  3549     merge_rls "simplify_rational" expand_binoms

  3550     (append_rls "divide" calculate_Rational

  3551 		[Thm ("real_divide_1",num_str real_divide_1),

  3552 		 (*"?x / 1 = ?x"*)

  3553 		 Thm ("rat_mult",num_str rat_mult),

  3554 		 (*(1)"?a / ?b * (?c / ?d) = ?a * ?c / (?b * ?d)"*)

  3555 		 Thm ("real_times_divide1_eq",num_str real_times_divide1_eq),

  3556 		 (*(2)"?a * (?c / ?d) = ?a * ?c / ?d" must be [2],

  3557 		 otherwise inv.to a / b / c = ...*)

  3558 		 Thm ("real_times_divide2_eq",num_str real_times_divide2_eq),

  3559 		 (*"?a / ?b * ?c = ?a * ?c / ?b"*)

  3560 		 Thm ("add_minus",num_str add_minus),

  3561 		 (*"?a + ?b - ?b = ?a"*)

  3562 		 Thm ("add_minus1",num_str add_minus1),

  3563 		 (*"?a - ?b + ?b = ?a"*)

  3564 		 Thm ("real_divide_minus1",num_str real_divide_minus1)

  3565 		 (*"?x / -1 = - ?x"*)

  3566 (*

  3567 ,

  3568 		 Thm ("",num_str )

  3569 *)

  3570 		 ]);

  3572 (*---------vvv-------------MG ab 1.07.2003--------------vvv-----------*)

  3574 (* ------------------------------------------------------------------ *)

  3575 (*                  Simplifier für beliebige Buchterme                *)

  3576 (* ------------------------------------------------------------------ *)

  3577 (*----------------------- norm_Rational_mg ---------------------------*)

  3578 (*. description of the simplifier see MG-DA.p.56ff .*)

  3579 (* ------------------------------------------------------------------- *)

  3580 val common_nominator_p_rls = prep_rls(

  3581   Rls {id = "common_nominator_p_rls", preconds = [],

  3582        rew_ord = ("dummy_ord",dummy_ord),

  3583 	 erls = e_rls, srls = Erls, calc = [],

  3584 	 rules =

  3585 	 [Rls_ common_nominator_p

  3586 	  (*FIXME.WN0401 ? redesign Rrls - use exhaustively on a term ?

  3587 	    FIXME.WN0510 unnecessary nesting: introduce RRls_ : rls -> rule*)

  3588 	  ],

  3589 	 scr = EmptyScr});

  3590 (* ------------------------------------------------------------------- *)

  3591 val cancel_p_rls = prep_rls(

  3592   Rls {id = "cancel_p_rls", preconds = [],

  3593        rew_ord = ("dummy_ord",dummy_ord),

  3594 	 erls = e_rls, srls = Erls, calc = [],

  3595 	 rules =

  3596 	 [Rls_ cancel_p

  3597 	  (*FIXME.WN.0401 ? redesign Rrls - use exhaustively on a term ?*)

  3598 	  ],

  3599 	 scr = EmptyScr});

  3600 (* -------------------------------------------------------------------- *)

  3601 (*. makes 'normal' fractions; 'is_polyexp' inhibits double fractions;

  3602     used in initial part norm_Rational_mg, see example DA-M02-main.p.60.*)

  3603 val rat_mult_poly = prep_rls(

  3604   Rls {id = "rat_mult_poly", preconds = [],

  3605        rew_ord = ("dummy_ord",dummy_ord),

  3606 	 erls =  append_rls "e_rls-is_polyexp" e_rls

  3607 	         [Calc ("Poly.is'_polyexp", eval_is_polyexp "")],

  3608 	 srls = Erls, calc = [],

  3609 	 rules =

  3610 	 [Thm ("rat_mult_poly_l",num_str rat_mult_poly_l),

  3611 	  (*"?c is_polyexp ==> ?c * (?a / ?b) = ?c * ?a / ?b"*)

  3612 	  Thm ("rat_mult_poly_r",num_str rat_mult_poly_r)

  3613 	  (*"?c is_polyexp ==> ?a / ?b * ?c = ?a * ?c / ?b"*)

  3614 	  ],

  3615 	 scr = EmptyScr});

  3616 (* ------------------------------------------------------------------ *)

  3617 (*. makes 'normal' fractions; 'is_polyexp' inhibits double fractions;

  3618     used in looping part norm_Rational_rls, see example DA-M02-main.p.60

  3619     .. WHERE THE LATTER DOES ALWAYS WORK, BECAUSE erls = e_rls,

  3620     I.E. THE RESPECTIVE ASSUMPTION IS STORED AND Thm APPLIED; WN051028

  3621     ... WN0609???MG.*)

  3622 val rat_mult_div_pow = prep_rls(

  3623   Rls {id = "rat_mult_div_pow", preconds = [],

  3624        rew_ord = ("dummy_ord",dummy_ord),

  3625        erls = e_rls,

  3626        (*FIXME.WN051028 append_rls "e_rls-is_polyexp" e_rls

  3627 			[Calc ("Poly.is'_polyexp", eval_is_polyexp "")],

  3628          with this correction ^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^ we get

  3629 	 error "rational.sml.sml: diff.behav. in norm_Rational_mg 29" etc.

  3630          thus we decided to go on with this flaw*)

  3631 		 srls = Erls, calc = [],

  3632       rules = [Thm ("rat_mult",num_str rat_mult),

  3633 	       (*"?a / ?b * (?c / ?d) = ?a * ?c / (?b * ?d)"*)

  3634 	       Thm ("rat_mult_poly_l",num_str rat_mult_poly_l),

  3635 	       (*"?c is_polyexp ==> ?c * (?a / ?b) = ?c * ?a / ?b"*)

  3636 	       Thm ("rat_mult_poly_r",num_str rat_mult_poly_r),

  3637 	       (*"?c is_polyexp ==> ?a / ?b * ?c = ?a * ?c / ?b"*)

  3639 	       Thm ("real_divide_divide1_mg", real_divide_divide1_mg),

  3640 	       (*"y ~= 0 ==> (u / v) / (y / z) = (u * z) / (y * v)"*)

  3641 	       Thm ("real_divide_divide1_eq", real_divide_divide1_eq),

  3642 	       (*"?x / (?y / ?z) = ?x * ?z / ?y"*)

  3643 	       Thm ("real_divide_divide2_eq", real_divide_divide2_eq),

  3644 	       (*"?x / ?y / ?z = ?x / (?y * ?z)"*)

  3645 	       Calc ("HOL.divide"  ,eval_cancel "#divide_"),

  3647 	       Thm ("rat_power", num_str rat_power)

  3648 		(*"(?a / ?b) ^^^ ?n = ?a ^^^ ?n / ?b ^^^ ?n"*)

  3649 	       ],

  3650       scr = Script ((term_of o the o (parse thy)) "empty_script")

  3651       }:rls);

  3652 (* ------------------------------------------------------------------ *)

  3653 val rat_reduce_1 = prep_rls(

  3654   Rls {id = "rat_reduce_1", preconds = [],

  3655        rew_ord = ("dummy_ord",dummy_ord),

  3656        erls = e_rls, srls = Erls, calc = [],

  3657        rules = [Thm ("real_divide_1",num_str real_divide_1),

  3658 		(*"?x / 1 = ?x"*)

  3659 		Thm ("real_mult_1",num_str real_mult_1)

  3660 		(*"1 * z = z"*)

  3661 		],

  3662        scr = Script ((term_of o the o (parse thy)) "empty_script")

  3663        }:rls);

  3664 (* ------------------------------------------------------------------ *)

  3665 (*. looping part of norm_Rational(*_mg*) .*)

  3666 val norm_Rational_rls = prep_rls(

  3667    Rls {id = "norm_Rational_rls", preconds = [],

  3668        rew_ord = ("dummy_ord",dummy_ord),

  3669        erls = norm_rat_erls, srls = Erls, calc = [],

  3670        rules = [Rls_ common_nominator_p_rls,

  3671 		Rls_ rat_mult_div_pow,

  3672 		Rls_ make_rat_poly_with_parentheses,

  3673 		Rls_ cancel_p_rls,(*FIXME:cancel_p does NOT order sometimes*)

  3674 		Rls_ rat_reduce_1

  3675 		],

  3676        scr = Script ((term_of o the o (parse thy)) "empty_script")

  3677        }:rls);

  3678 (* ------------------------------------------------------------------ *)

  3679 (*040109 'norm_Rational'(by RL) replaced by 'norm_Rational_mg'(MG)

  3680  just be renaming:*)

  3681 val norm_Rational(*_mg*) = prep_rls(

  3682    Seq {id = "norm_Rational"(*_mg*), preconds = [],

  3683        rew_ord = ("dummy_ord",dummy_ord),

  3684        erls = norm_rat_erls, srls = Erls, calc = [],

  3685        rules = [Rls_ discard_minus_,

  3686 		Rls_ rat_mult_poly,(* removes double fractions like a/b/c    *)

  3687 		Rls_ make_rat_poly_with_parentheses, (*WN0510 also in(#)below*)

  3688 		Rls_ cancel_p_rls, (*FIXME.MG:cancel_p does NOT order sometim*)

  3689 		Rls_ norm_Rational_rls,   (* the main rls, looping (#)       *)

  3690 		Rls_ discard_parentheses_ (* mult only                       *)

  3691 		],

  3692        scr = Script ((term_of o the o (parse thy)) "empty_script")

  3693        }:rls);

  3694 (* ------------------------------------------------------------------ *)

  3697 ruleset' := overwritelthy thy (!ruleset',

  3698   [("calculate_Rational", calculate_Rational),

  3699    ("calc_rat_erls",calc_rat_erls),

  3700    ("rational_erls", rational_erls),

  3701    ("cancel_p", cancel_p),

  3702    ("cancel", cancel),

  3703    ("common_nominator_p", common_nominator_p),

  3704    ("common_nominator_p_rls", common_nominator_p_rls),

  3705    ("common_nominator"  , common_nominator),

  3706    ("discard_minus", discard_minus),

  3707    ("powers_erls", powers_erls),

  3708    ("powers", powers),

  3709    ("rat_mult_divide", rat_mult_divide),

  3710    ("reduce_0_1_2", reduce_0_1_2),

  3711    ("rat_reduce_1", rat_reduce_1),

  3712    ("norm_rat_erls", norm_rat_erls),

  3713    ("norm_Rational", norm_Rational),

  3714    ("norm_Rational_rls", norm_Rational_rls),

  3715    ("norm_Rational_parenthesized", norm_Rational_parenthesized),

  3716    ("rat_mult_poly", rat_mult_poly),

  3717    ("rat_mult_div_pow", rat_mult_div_pow),

  3718    ("cancel_p_rls", cancel_p_rls)

  3719    ]);

  3721 calclist':= overwritel (!calclist',

  3722    [("is_expanded", ("Rational.is'_expanded", eval_is_expanded ""))

  3723     ]);

  3725 (** problems **)

  3727 store_pbt

  3728  (prep_pbt Rational.thy "pbl_simp_rat" [] e_pblID

  3729  (["rational","simplification"],

  3730   [("#Given" ,["term t_"]),

  3731    ("#Where" ,["t_ is_ratpolyexp"]),

  3732    ("#Find"  ,["normalform n_"])

  3733   ],

  3734   append_rls "e_rls" e_rls [(*for preds in where_*)],

  3735   Some "Simplify t_",

  3736   [["simplification","of_rationals"]]));

  3738 (** methods **)

  3740 store_met

  3741     (prep_met Rational.thy "met_simp_rat" [] e_metID

  3742 	      (["simplification","of_rationals"],

  3743 	       [("#Given" ,["term t_"]),

  3744 		("#Where" ,["t_ is_ratpolyexp"]),

  3745 		("#Find"  ,["normalform n_"])

  3746 		],

  3747 	       {rew_ord'="tless_true",

  3748 		rls' = e_rls,

  3749 		calc = [], srls = e_rls,

  3750 		prls = append_rls "simplification_of_rationals_prls" e_rls

  3751 				[(*for preds in where_*)

  3752 				 Calc ("Rational.is'_ratpolyexp",

  3753 				       eval_is_ratpolyexp "")],

  3754 		crls = e_rls, nrls = e_rls},

  3755 	       "Script SimplifyScript (t_::real) =                \

  3756 	       \  ((Rewrite_Set norm_Rational False) t_)"

  3757 	       ));

  3759 (* use"../IsacKnowledge/Rational.ML";

  3760    use"IsacKnowledge/Rational.ML";

  3761    use"Rational.ML";

  3762    *)

author	wneuper
	Thu, 14 Sep 2006 15:29:45 +0200
branch	start_Take
changeset 660	c58b3d2ee7a9
parent 650	4a9ff889afaf
child 669	aaa641d8546d
permissions	-rw-r--r--