wneuper/isa: src/Tools/isac/Knowledge/Rational.thy@02a9b2540eb7

     1 (* rationals, i.e. fractions of multivariate polynomials over the real field

     2    author: isac team

     3    Copyright (c) isac team 2002

     4    Use is subject to license terms.

     6    depends on Poly (and not on Atools), because

     7    fractions with _normalized_ polynomials are canceled, added, etc.

     8 *)

    10 theory Rational imports Poly begin

    12 consts

    14   is'_expanded   :: "real => bool" ("_ is'_expanded")     (*RL->Poly.thy*)

    15   is'_ratpolyexp :: "real => bool" ("_ is'_ratpolyexp")

    17 axioms (*.not contained in Isabelle2002,

    18           stated as axioms, TODO?: prove as theorems*)

    20   mult_cross:      "[| b ~= 0; d ~= 0 |] ==> (a / b = c / d) = (a * d = b * c)"

    21   mult_cross1:     "   b ~= 0            ==> (a / b = c    ) = (a     = b * c)"

    22   mult_cross2:     "           d ~= 0    ==> (a     = c / d) = (a * d =     c)"

    24   add_minus:       "a + b - b = a"(*RL->Poly.thy*)

    25   add_minus1:      "a - b + b = a"(*RL->Poly.thy*)

    27   rat_mult:        "a / b * (c / d) = a * c / (b * d)"(*?Isa02*)

    28   rat_mult2:       "a / b *  c      = a * c /  b     "(*?Isa02*)

    30   rat_mult_poly_l: "c is_polyexp ==> c * (a / b) = c * a /  b"

    31   rat_mult_poly_r: "c is_polyexp ==> (a / b) * c = a * c /  b"

    33 (*real_times_divide1_eq .. Isa02*)

    34   real_times_divide_1_eq:  "-1    * (c / d) =-1 * c /      d "

    35   real_times_divide_num:   "a is_const ==>

    36 	          	   a     * (c / d) = a * c /      d "

    38   real_mult_div_cancel2:   "k ~= 0 ==> m * k / (n * k) = m / n"

    39 (*real_mult_div_cancel1:   "k ~= 0 ==> k * m / (k * n) = m / n"..Isa02*)

    41   real_divide_divide1:     "y ~= 0 ==> (u / v) / (y / z) = (u / v) * (z / y)"

    42   real_divide_divide1_mg:  "y ~= 0 ==> (u / v) / (y / z) = (u * z) / (y * v)"

    43 (*real_divide_divide2_eq:  "x / y / z = x / (y * z)"..Isa02*)

    45   rat_power:               "(a / b)^^^n = (a^^^n) / (b^^^n)"

    48   rat_add:         "[| a is_const; b is_const; c is_const; d is_const |] ==>

    49 	           a / c + b / d = (a * d + b * c) / (c * d)"

    50   rat_add_assoc:   "[| a is_const; b is_const; c is_const; d is_const |] ==>

    51 	           a / c +(b / d + e) = (a * d + b * c)/(d * c) + e"

    52   rat_add1:        "[| a is_const; b is_const; c is_const |] ==>

    53 	           a / c + b / c = (a + b) / c"

    54   rat_add1_assoc:   "[| a is_const; b is_const; c is_const |] ==>

    55 	           a / c + (b / c + e) = (a + b) / c + e"

    56   rat_add2:        "[| a is_const; b is_const; c is_const |] ==>

    57 	           a / c + b = (a + b * c) / c"

    58   rat_add2_assoc:  "[| a is_const; b is_const; c is_const |] ==>

    59 	           a / c + (b + e) = (a + b * c) / c + e"

    60   rat_add3:        "[| a is_const; b is_const; c is_const |] ==>

    61 	           a + b / c = (a * c + b) / c"

    62   rat_add3_assoc:   "[| a is_const; b is_const; c is_const |] ==>

    63 	           a + (b / c + e) = (a * c + b) / c + e"

    65 text {*calculate in rationals: gcd, lcm, etc.

    66       (c) Stefan Karnel 2002

    67       Institute for Mathematics D and Institute for Software Technology,

    68       TU-Graz SS 2002 *}

    70 text {* Remark on notions in the documentation below:

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

    72     [2] 'polynomial' normalform (Polynom)

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

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

    75     Instead of

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

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

    78     we say:

    79       'fraction in normalform [2]'

    80       'fraction in normalform [3]'

    81     or

    82       'fraction [2]'

    83       'fraction [3]'.

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

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

    86 *}

    88 ML {*

    89 val thy = @{theory};

    91 signature RATIONALI =

    92 sig

    93   type mv_monom

    94   type mv_poly

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

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

    97   val calculate_Rational : rls

    98   val calc_rat_erls:rls

    99   val cancel : rls

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

   101   val cancel_p : rls

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

   103   val common_nominator : rls

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

   105   val common_nominator_p : rls

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

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

   108 			 (string * term) option

   109   val expanded2polynomial : term -> term option

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

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

   112   val is_expanded : term -> bool

   113   val is_polynomial : term -> bool

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

   116   val mv_lcm : mv_poly -> mv_poly -> mv_poly

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

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

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

   121                            roots into a simple and canceled fraction

   122                            with normalform [2].*)

   123 *)

   124 (*val norm_rational_p : 19.10.02 missing FIXXXXXXXXXXXXME

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

   126                            roots into a simple and canceled fraction

   127                            with normalform [2].*)

   128 *)

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

   130   val polynomial2expanded : term -> term option

   131   val rational_erls :

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

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

   135 end

   137 (*.**************************************************************************

   138 survey on the functions

   139 ~~~~~~~~~~~~~~~~~~~~~~~

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

   141 --------------------:------------------+-------------------:-----------------

   142  factout_p_         :                  | factout_          :

   143  cancel_p_          :                  | cancel_           :

   144                     :cancel_p          |                   :cancel

   145 --------------------:------------------+-------------------:-----------------

   146  common_nominator_p_:                  | common_nominator_ :

   147                     :common_nominator_p|                   :common_nominator

   148  add_fraction_p_    :                  | add_fraction_     :

   149 --------------------:------------------+-------------------:-----------------

   150 ???SK                 :norm_rational_p   |                   :norm_rational

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

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

   156 conversion from Isabelle-term to internal representation

   157 ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~

   159 ... BITTE FORTSETZEN ...

   161 polynomial2expanded = ...

   162 expanded2polynomial = ...

   164 remark: polynomial2expanded o expanded2polynomial = I,

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

   167 functions for greatest common divisor and canceling

   168 ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~

   169 ################################################################################

   170 ##   search Isabelle2009-2/src/HOL/Multivariate_Analysis

   171 ##   Amine Chaieb, Robert Himmelmann, John Harrison.

   172 ################################################################################

   173 mv_gcd

   174 factout_

   175 factout_p_

   176 cancel_

   177 cancel_p_

   179 functions for least common multiple and addition of fractions

   180 ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~

   181 mv_lcm

   182 common_nominator_

   183 common_nominator_p_

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

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

   187 functions for normalform of rationals

   188 ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~

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

   190           schaffen diese SML-Funktionen wirklich ganz allgemeine Terme ?

   192 norm_rational_

   193 norm_expanded_rat_

   195 **************************************************************************.*)

   198 (*##*)

   199 structure RationalI : RATIONALI =

   200 struct

   201 (*##*)

   203 infix mem ins union; (*WN100819 updating to Isabelle2009-2*)

   204 fun x mem [] = false

   205   | x mem (y :: ys) = x = y orelse x mem ys;

   206 fun (x ins xs) = if x mem xs then xs else x :: xs;

   207 fun xs union [] = xs

   208   | [] union ys = ys

   209   | (x :: xs) union ys = xs union (x ins ys);

   211 (*. gcd of integers .*)

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

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

   214 		  else gcd_int b (a mod b);

   216 (*. univariate polynomials (uv) .*)

   217 (*. univariate polynomials are represented as a list

   218     of the coefficent in reverse maximum degree order .*)

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

   220 type uv_poly = int list;

   222 (*. adds two uv polynomials .*)

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

   224   | uv_mod_add_poly (p1,[]) = p1

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

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

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

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

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

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

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

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

   236 fun uv_mod_deg ([]:uv_poly) = 0

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

   239 (*. calculates the remainder of x/p and represents it as

   240     value between -p/2 and p/2 .*)

   241 fun uv_mod_mod2(x,p)=

   242     let

   243 	val y=(x mod p);

   244     in

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

   246 	    (

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

   248 	     )

   249     end;

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

   252 fun uv_mod_list_modp [] p = []

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

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

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

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

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

   260 (*. only for uv_mod_divides .*)

   261 (* FIXME: Division von x^9+x^5+1 durch x-1000 funktioniert nicht,

   262    integer zu klein  *)

   263 fun uv_mod_pdiv (p1:uv_poly) ([]:uv_poly) =

   264     error ("RATIONALS_UV_MOD_PDIV_EXCEPTION: division by zero")

   265   | uv_mod_pdiv p1 [x] =

   266     let

   267 	val xs= Unsynchronized.ref  [];

   268     in

   269 	if x<>0 then

   270 	    (

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

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

   273 	     )

   274 	else error ("RATIONALS_UV_MOD_PDIV_EXCEPTION: division by zero");

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

   276     end

   277   | uv_mod_pdiv p1 p2 =

   278     let

   279 	val n= uv_mod_deg(p2);

   280 	val m= Unsynchronized.ref (uv_mod_deg(p1));

   281 	val p1'= Unsynchronized.ref  (rev(p1));

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

   283 	val lc2=hd(p2');

   284 	val q= Unsynchronized.ref  [];

   285 	val c= Unsynchronized.ref  0;

   286 	val output= Unsynchronized.ref  ([],[]);

   287     in

   288 	(

   289 	 if (!m)=0 orelse p2=[0]

   290          then error ("RATIONALS_UV_MOD_PDIV_EXCEPTION: Division by zero")

   291 	 else

   292 	     (

   293 	      if (!m)<n then

   294 		  (

   295 		   output:=([0],p1)

   296 		   )

   297 	      else

   298 		  (

   299 		   while (!m)>=n do

   300 		       (

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

   302 			if !c<>0 then

   303 			    (

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

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

   306 			     m:=uv_mod_deg(!p1')

   307 			     )

   308 			else m:=0

   309 			);

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

   311 		   )

   312 	      );

   313 	     !output

   314 	 )

   315     end;

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

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

   319     let

   320 	val n=uv_mod_deg(p2);

   321 	val m= Unsynchronized.ref  (uv_mod_deg(uv_mod_list_modp p1 p));

   322 	val p1'= Unsynchronized.ref  (rev(p1));

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

   324 	val lc2=hd(p2');

   325 	val q= Unsynchronized.ref  [];

   326 	val c= Unsynchronized.ref  0;

   327 	val output= Unsynchronized.ref  ([],[]);

   328     in

   329 	(

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

   331 	 else

   332 	     (

   333 	      if (!m)<n then

   334 		  (

   335 		   output:=([0],p1)

   336 		   )

   337 	      else

   338 		  (

   339 		   while (!m)>=n do

   340 		       (

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

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

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

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

   345 			m:=(!m)-1

   346 			);

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

   349 		       (

   350 			p1':=tl(!p1')

   351 			);

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

   354 		   )

   355 	      );

   356 	     !output:uv_poly * uv_poly

   357 	 )

   358     end;

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

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

   362   | uv_mod_prest [] p2 = []:uv_poly

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

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

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

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

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

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

   371 fun uv_mod_cont ([]:uv_poly) = 0

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

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

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

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

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

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

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

   381   | uv_mod_pp p =

   382     let

   383 	val c= Unsynchronized.ref  0;

   384     in

   385 	(

   386 	 c:=uv_mod_cont(p);

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

   389 	 else uv_mod_div_list(p,!c)

   390 	)

   391     end;

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

   394 fun uv_mod_lc ([]:uv_poly) = 0

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

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

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

   399     let

   400 	val f = Unsynchronized.ref  [];

   401 	val f'= Unsynchronized.ref  p2;

   402 	val fi= Unsynchronized.ref  [];

   403     in

   404 	(

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

   406  	 while uv_mod_deg(!f')>0 do

   407 	     (

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

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

   410 		  (

   411 		   fi:=(!f');

   412 		   f:=(!fi)::(!f)

   413 		   )

   414 	      else ()

   415 	      );

   416 	      (!f)

   418 	 )

   419     end;

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

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

   423   | uv_mod_gcd_modp p1 [] p= p1

   424   | uv_mod_gcd_modp p1 p2 p=

   425     let

   426 	val p1'= Unsynchronized.ref [];

   427 	val p2'= Unsynchronized.ref [];

   428 	val pc= Unsynchronized.ref [];

   429 	val g= Unsynchronized.ref  [];

   430 	val d= Unsynchronized.ref  0;

   431 	val prs= Unsynchronized.ref  [];

   432     in

   433 	(

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

   435 	     (

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

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

   438 	      )

   439 	 else

   440 	     (

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

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

   443 	      );

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

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

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

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

   450 	 else pc:=hd(!prs);

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

   453 	 !g

   454 	 )

   455     end;

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

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

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

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

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

   464 fun uv_mod_add_qu [] = 0.0

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

   467 (*. calculates the euklidean norm .*)

   468 fun uv_mod_norm ([]:uv_poly) = 0.0

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

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

   472 fun uv_mod_multi a b = a * b;

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

   475 fun uv_mod_prim(x,[])= false

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

   477 		else false

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

   479 			then

   480 			    if uv_mod_prim(x,ys) then true

   481 			    else false

   482 		    else false;

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

   485 fun uv_mod_nextprime(g,p)=

   486     let

   487 	val list= Unsynchronized.ref  [2];

   488 	val exit= Unsynchronized.ref  0;

   489 	val i = Unsynchronized.ref 2

   490     in

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

   492 	    (

   493 	     if uv_mod_prim(!i,!list) then

   494 		 (

   495 		  if (p mod !i <> 0)

   496 		      then

   497 			  (

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

   499 			   i:= (!i)+1

   500 			   )

   501 		  else i:=(!i)+1

   502 		  )

   503 	     else i:= (!i)+1

   504 		 );

   505 	    i:=(p+1);

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

   507 	    (

   508 	     if uv_mod_prim(!i,!list) then

   509 		 (

   510 		  if (g mod !i <> 0)

   511 		      then

   512 			  (

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

   514 			   exit:= (!i)

   515 			   )

   516 		  else i:=(!i)+1

   517 		      )

   518 	     else i:= (!i)+1

   519 		 );

   520 	    !exit

   521     end;

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

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

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

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

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

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

   531 (*. chinese remainder algorithm .*)

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

   533     let

   534 	val c= Unsynchronized.ref  0;

   535 	val r1'= Unsynchronized.ref  0;

   536 	val d= Unsynchronized.ref  0;

   537 	val a= Unsynchronized.ref  0;

   538     in

   539 	(

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

   541 	     (

   542 	      c:=(!c)+1

   543 	      );

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

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

   546 	 !r1'+(!d)*m1

   547 	 )

   548     end;

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

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

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

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

   554   | 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));

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

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

   558     let

   559 	val p1= Unsynchronized.ref  (uv_mod_pp(p1'));

   560 	val p2= Unsynchronized.ref  (uv_mod_pp(p2'));

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

   562 	val temp= Unsynchronized.ref  [];

   563 	val cp= Unsynchronized.ref  [];

   564 	val qp= Unsynchronized.ref  [];

   565 	val q= Unsynchronized.ref [];

   566 	val pn= Unsynchronized.ref  0;

   567 	val d= Unsynchronized.ref  0;

   568 	val g1= Unsynchronized.ref  0;

   569 	val p= Unsynchronized.ref  0;

   570 	val m= Unsynchronized.ref  0;

   571 	val exit= Unsynchronized.ref  0;

   572 	val i= Unsynchronized.ref  1;

   573     in

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

   575 	else

   576 	    (

   577 	     temp:= !p1;

   578 	     p1:= !p2;

   579 	     p2:= !temp

   580 	     );

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

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

   585 	p:=4;

   587 	m:=Real.ceil(2.0 * Real.fromInt(!d) *

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

   589 	  Real.fromInt(!d) *

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

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

   593 	while (!exit=0) do

   594 	    (

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

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

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

   598 		 (

   599 		  i:=1;

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

   601 		      (

   602 		       i:=(!i)+1

   603 		       );

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

   605 		  )

   606 	     else ();

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

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

   612 	     pn:=(!p);

   613 	     q:=(!qp);

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

   616 		 (

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

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

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

   620  		      (

   621  		       i:=1;

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

   623  			   (

   624  			    i:=(!i)+1

   625 		           );

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

   627  		      )

   628  		  else ();

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

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

   632  		      (

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

   634  		       )

   635  		  else

   636  		      (

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

   638  			   (

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

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

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

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

   643 		 	    )

   644 		       else

   645 			   (

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

   647 				(

   648 				 pn:= !p;

   649 				 q:= !qp

   650 				 )

   651 			    else ()

   652 			    )

   653 		       )

   654 		  );

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

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

   657 	     );

   658 	    uv_mod_smul_poly(!q,c):uv_poly

   659     end;

   661 (*. multivariate polynomials .*)

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

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

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

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

   667 (*. global variables .*)

   668 (*. order indicators .*)

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

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

   672 (*. datatypes for internal representation.*)

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

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

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

   677 type mv_poly = mv_monom list;

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

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

   681 fun mv_mg_hlp([]) = EQUAL

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

   683 		    else if x>0 then GREATER

   684 			 else mv_mg_hlp(list);

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

   687 fun mv_addlist([]) = 0

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

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

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

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

   693     if order=LEX_ then

   694 	(

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

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

   697 	     )

   698     else

   699 	if order=GGO_ then

   700 	    (

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

   702 	     else

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

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

   705 	     )

   706 	else error ("RATIONALS_MV_MONOM_GREATER_EXCEPTION: Wrong Order");

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

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

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

   711 let

   712     val temp= Unsynchronized.ref  EQUAL;

   713 in

   714     if order=LEX_ then

   715 	(

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

   717 	     error ("RATIONALS_MV_GEQ_EXCEPTION: Order error")

   718 	 else

   719 	     (

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

   721 	      if !temp=EQUAL then

   722 		  ( if x1=x2 then EQUAL

   723 		    else if x1>x2 then GREATER

   724 			 else LESS

   725 			     )

   726 	      else (!temp)

   727 	      )

   728 	     )

   729     else

   730 	if order=GGO_ then

   731 	    (

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

   733 	      error ("RATIONALS_MV_GEQ_EXCEPTION: Order error")

   734 	     else

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

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

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

   738 		     )

   739 	else error ("RATIONALS_MV_GEQ_EXCEPTION: Wrong Order")

   740 end;

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

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

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

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

   747 (*. leading power product .*)

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

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

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

   752 (*. leading monomial .*)

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

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

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

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

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

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

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

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

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

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

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

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

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

   770   | mv_lc(p1,order)  =

   771     let

   772 	val p1o= Unsynchronized.ref  (rev(sort (mv_geq order) (mv_rev_to(p1))));

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

   774 	val lc= Unsynchronized.ref  [];

   775     in

   776 	(

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

   778 	     (

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

   780 	      p1o:=tl(!p1o)

   781 	      );

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

   783 	 mv_rev_to(!lc)

   784 	 )

   785     end;

   787 (*. compares two powerproducts .*)

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

   790 (*. help function for mv_add .*)

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

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

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

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

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

   796   | mv_madd(p1,p2,order) =

   797     (

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

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

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

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

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

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

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

   805 	      )

   807 (*. adds two multivariate polynomials .*)

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

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

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

   812 (*. monom multiplication .*)

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

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

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

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

   819 (*. zeros a list .*)

   820 fun mv_null2([])=[]

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

   823 (*. multiplies two multivariate polynomials .*)

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

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

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

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

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

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

   831 fun mv_getmax([])=0

   832   | mv_getmax(x::p1)= let

   833 		       val m=mv_getmax(p1);

   834 		   in

   835 		       if m>x then m

   836 		       else x

   837 		   end;

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

   839 fun mv_grad([]:mv_poly) = 0

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

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

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

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

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

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

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

   850 fun mv_is_negativ([])=false

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

   853 (*. division of monomials .*)

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

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

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

   857     let

   858 	val c= Unsynchronized.ref  (#1(p2));

   859 	val pp= Unsynchronized.ref  [];

   860     in

   861 	(

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

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

   864 	     if #1(p2)<>0 then

   865 		 (

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

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

   868 		  else (!c,!pp)

   869 		      )

   870 	     else error("MV_MDIV_EXCEPTION Dividing by empty Polynom")

   871 		 )

   872     end;

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

   875 fun mv_print_poly([]:mv_poly)=tracing("[]\n")

   876   | mv_print_poly((x,y)::[])= tracing("("^Int.toString(x)^","^ints2str(y)^")\n")

   877   | mv_print_poly((x,y)::p1) = (tracing("("^Int.toString(x)^","^ints2str(y)^"),");mv_print_poly(p1));

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

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

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

   883   | mv_division(f,g,order)=

   884     let

   885 	val r= Unsynchronized.ref  [];

   886 	val q= Unsynchronized.ref  [];

   887 	val g'= Unsynchronized.ref  ([] : mv_monom list);

   888 	val k= Unsynchronized.ref  0;

   889 	val m= Unsynchronized.ref  (0,[0]);

   890 	val exit= Unsynchronized.ref  0;

   891     in

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

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

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

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

   896 	else

   897 	    (

   898 	     exit:=0;

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

   900 		 (

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

   902 		  else error ("RATIONALS_MV_DIVISION_EXCEPTION: Dividing by zero");

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

   904 		      (

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

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

   907 		       )

   908 		  else exit:=1;

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

   910 		  else (exit:=1)

   911 		  );

   912 		 (rev(!q),!r)

   913 		 )

   914     end;

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

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

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

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

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

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

   924 (*. multivariate case .*)

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

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

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

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

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

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

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

   935 fun mv_deg([]:mv_poly) = 0

   936   | mv_deg(p1)=

   937     let

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

   939     in

   940 	if length(p1')=0 then 0

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

   942     end;

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

   945 fun mv_deg2([]:mv_poly) = 0

   946   | mv_deg2(p1)=

   947     let

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

   949     in

   950 	if length(p1')=0 then 0

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

   952     end;

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

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

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

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

   959 fun uv_content2([]:mv_poly) = 0

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

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

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

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

   965     let

   966 	val count= Unsynchronized.ref  0;

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

   968 	val help= Unsynchronized.ref  ((c1,e1)::others);

   969 	val list= Unsynchronized.ref  [];

   970     in

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

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

   973 	     else

   974 		 (

   975 		  count:=0;

   976 		  while (!count)<=max do

   977 		      (

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

   979 			   (

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

   981 			    help:=tl(!help)

   982 			    )

   983 		       else

   984 			   (

   985 			    list:= 0::(!list)

   986 			    );

   987 		       count := (!count) + 1

   988 		       );

   989 		      (!list)

   990 		      )

   991     end;

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

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

   995   | uv_to_poly p1 =

   996     let

   997 	val count= Unsynchronized.ref  0;

   998 	val help= Unsynchronized.ref  p1;

   999 	val list= Unsynchronized.ref  [];

  1000     in

  1001 	while length(!help)>0 do

  1002 	    (

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

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

  1005 	     count:=(!count)+1;

  1006 	     help:=tl(!help)

  1007 	     );

  1008 	    (!list)

  1009     end;

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

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

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

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

  1015     let

  1016 	val list= Unsynchronized.ref  (rev(sort (mv_geq LEX_) (mv_shorten(p1,LEX_))));

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

  1018     in

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

  1020     end

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

  1022     let

  1023 	val list= Unsynchronized.ref  (rev(sort (mv_geq LEX_) (mv_shorten(p2,LEX_))));

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

  1025     in

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

  1027     end

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

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

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

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

  1033     let

  1034 	val x= Unsynchronized.ref  (rev(pl));

  1035 	val t= Unsynchronized.ref  [];

  1036 	val y= Unsynchronized.ref  [];

  1037 	val n= Unsynchronized.ref  1;

  1038 	val n1= Unsynchronized.ref [];

  1039     in

  1040 	(

  1041 	 while length(!x)>1 do

  1042 	     (

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

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

  1045 		   else n1:=[];

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

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

  1048 	      x:=tl(!x)

  1049 	      );

  1050 	 (!y)

  1051 	 )

  1052     end;

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

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

  1056   | mv_newton_add [x:mv_poly] t = x

  1057   | mv_newton_add (pl:mv_poly list) t =

  1058     let

  1059 	val expos= Unsynchronized.ref  [];

  1060 	val pll= Unsynchronized.ref  pl;

  1061     in

  1062 	(

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

  1065 	     (

  1066 	      pll:=tl(!pll)

  1067 	      );

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

  1069 	 mv_add(hd(pl),

  1070 		mv_mul(

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

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

  1073 		       LEX_

  1074 		       ),

  1075 		LEX_)

  1076 	 )

  1077     end;

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

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

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

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

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

  1084   | mv_newton pl =

  1085     let

  1086 	val c= Unsynchronized.ref  pl;

  1087 	val c1= Unsynchronized.ref  [];

  1088 	val n=length(pl);

  1089 	val k= Unsynchronized.ref  1;

  1090 	val i= Unsynchronized.ref  n;

  1091 	val ppl= Unsynchronized.ref  [];

  1092     in

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

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

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

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

  1097 	    (

  1098 	     k:=(!k)+1;

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

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

  1101 	     ppl:= !c;

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

  1103 	     );

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

  1105 	c1:=rev(!c1);

  1106 	ppl:= !c1;

  1107 	mv_newton_add (!c1) 1

  1108     end;

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

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

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

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

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

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

  1117     let

  1118 	val y= Unsynchronized.ref  xs;

  1119 	val x= Unsynchronized.ref  [];

  1120     in

  1121 	(

  1122 	 x:=e;

  1123 	 while length(!y)>0 do

  1124 	     (

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

  1126 	      y:=tl(!y)

  1127 	      );

  1128 	 !x

  1129 	 )

  1130     end;

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

  1133 fun list_is_null [] = true

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

  1136 (* check if main variable is zero*)

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

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

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

  1141   | mv_content(p1) =

  1142     let

  1143 	val list= Unsynchronized.ref  (rev(sort (mv_geq LEX_) (mv_shorten(p1,LEX_))));

  1144 	val test= Unsynchronized.ref  (hd(#2(hd(!list))));

  1145 	val result= Unsynchronized.ref  [];

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

  1147     in

  1148 	(

  1149 	 if length(!list)>1 then

  1150 	     (

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

  1152 		  (

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

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

  1156 		   else list:=tl(!list)

  1158 		       );

  1159 		  if length(!list)>0 then

  1160 		   (

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

  1162 		    )

  1163 		  else list:=(!result);

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

  1165 		  (!list)

  1166 		  )

  1167 	 else

  1168 	     (

  1169 	      mv_null3(!list)

  1170 	      )

  1171 	     )

  1172     end

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

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

  1176   | mv_pp(p1) = let

  1177 		    val cont= Unsynchronized.ref  [];

  1178 		    val pp= Unsynchronized.ref [];

  1179 		in

  1180 		    cont:=mv_content(p1);

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

  1182 		    if !pp=[]

  1183 			then error("RATIONALS_MV_PP_EXCEPTION: Invalid Content ")

  1184 		    else (!pp)

  1185 		end

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

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

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

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

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

  1192      let

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

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

  1195      in

  1196 	(

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

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

  1199         )

  1200     end

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

  1202 	(

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

  1204 	)

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

  1206 	(

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

  1208         )

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

  1210     let

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

  1212 	val cont =

  1213 		  (

  1214                    if main_zero(mv_content(p1')) andalso

  1215                      (main_zero(mv_content(p2'))) then

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

  1217                    else

  1218                      mv_gcd (mv_content(p1')) (mv_content(p2'))

  1219                   );

  1220 	val p1= #1(mv_division(p1',mv_content(p1'),LEX_));

  1221 	val p2= #1(mv_division(p2',mv_content(p2'),LEX_));

  1222 	val gcd= Unsynchronized.ref  [];

  1223 	val candidate= Unsynchronized.ref  [];

  1224 	val interpolation_list= Unsynchronized.ref  [];

  1225 	val delta= Unsynchronized.ref  [];

  1226         val p1r = Unsynchronized.ref [];

  1227         val p2r = Unsynchronized.ref [];

  1228         val p1r' = Unsynchronized.ref [];

  1229         val p2r' = Unsynchronized.ref [];

  1230 	val factor= Unsynchronized.ref  [];

  1231 	val r= Unsynchronized.ref  0;

  1232 	val gcd_r= Unsynchronized.ref  [];

  1233 	val d= Unsynchronized.ref  0;

  1234 	val exit= Unsynchronized.ref  0;

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

  1236     in

  1237 	(

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

  1239 	     (

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

  1241 	      )

  1242 	 else

  1243 	     (

  1244 	      while !exit=0 do

  1245 		  (

  1246 		   r:=(!r)+1;

  1247                    p1r := mv_lc(p1,LEX_);

  1248 		   p2r := mv_lc(p2,LEX_);

  1249                    if main_zero(!p1r) andalso

  1250                       main_zero(!p2r)

  1251                    then

  1252                        (

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

  1254                        )

  1255                    else

  1256                        (

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

  1258                        );

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

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

  1261 		   if mv_lc2(mv_shorten(mv_subs(!p1r,!r),LEX_),LEX_)=0 andalso

  1262 		      mv_lc2(mv_shorten(mv_subs(!p2r,!r),LEX_),LEX_)=0

  1263                    then

  1264                        (

  1265 		       )

  1266 		   else

  1267 		       (

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

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

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

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

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

  1273 			if (!d < !current_degree) then

  1274 			    (

  1275 			     current_degree:= !d;

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

  1277 			     )

  1278 			else

  1279 			    (

  1280 			     if (!d = !current_degree) then

  1281 				 (

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

  1283 				  )

  1284 			     else ()

  1285 				 )

  1286 			    );

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

  1288 			  (

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

  1290 			   if !candidate=[] then ()

  1291 			   else

  1292 			       (

  1293 				candidate:=mv_pp(!candidate);

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

  1295 				    (

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

  1297 				     exit:=1

  1298 				     )

  1299 				else ()

  1300 				    );

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

  1302 			       )

  1303 		      else ()

  1304 			  )

  1305 	     );

  1306 	     (!gcd):mv_poly

  1307 	     )

  1308     end;

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

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

  1313     (

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

  1315      );

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

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

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

  1321 fun count_neg ([]:mv_poly) = 0

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

  1323 			  else count_neg xs;

  1325 (*. help function for is_polynomial

  1326     checks the order of the operators .*)

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

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

  1329   | test_polynomial (t as Const ("Groups.times_class.times",_) $ t1 $ t2) v = if v="^" then false

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

  1331   | test_polynomial (t as Const ("Groups.plus_class.plus",_) $ t1 $ t2) v = if v="*" orelse v="^" then false

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

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

  1334   | test_polynomial _ v = false;

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

  1337 fun is_polynomial t = test_polynomial t " ";

  1339 (*. help function for is_expanded

  1340     checks the order of the operators .*)

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

  1342   | test_exp (t as Const ("Groups.times_class.times",_) $ t1 $ t2) v = if v="^" then false

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

  1344   | test_exp (t as Const ("Groups.plus_class.plus",_) $ t1 $ t2) v = if v="*" orelse v="^" then false

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

  1346   | test_exp (t as Const ("Groups.minus_class.minus",_) $ t1 $ t2) v = if v="*" orelse v="^" then false

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

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

  1349   | test_exp  _ v = false;

  1352 (*. help function for check_coeff:

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

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

  1355     let

  1356 	val x= Unsynchronized.ref  NONE;

  1357 	val len= Unsynchronized.ref  0;

  1358 	val vl= Unsynchronized.ref  [];

  1359 	val vh= Unsynchronized.ref  [];

  1360 	val i= Unsynchronized.ref  0;

  1361     in

  1362 	if is_numeral str then

  1363 	    (

  1364 	     SOME [(((the o int_of_str) str),mv_null2(v))] handle _ => NONE

  1365 		 )

  1366 	else (* variable *)

  1367 	    (

  1368 	     len:=length(v);

  1369 	     vh:=v;

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

  1371 		 (

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

  1373 		      (

  1374 		       vl:=1::(!vl)

  1375 		       )

  1376 		  else

  1377 		      (

  1378 		       vl:=0::(!vl)

  1379 		       );

  1380 		      vh:=tl(!vh);

  1381 		      i:=(!i)+1

  1382 		      );

  1383 		 SOME [(1,rev(!vl))] handle _ => NONE

  1384 	    )

  1385     end

  1386   | term2coef' (Const ("Groups.times_class.times",_) $ t1 $ t2) v :mv_poly option=

  1387     let

  1388 	val t1pp= Unsynchronized.ref  [];

  1389 	val t2pp= Unsynchronized.ref  [];

  1390 	val t1c= Unsynchronized.ref  0;

  1391 	val t2c= Unsynchronized.ref  0;

  1392     in

  1393 	(

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

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

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

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

  1399 	 SOME [( (!t1c)*(!t2c) ,( (map op+) ((!t1pp)~~(!t2pp)) ) )] handle _ => NONE

  1401 	 )

  1402     end

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

  1404     let

  1405 	val x= Unsynchronized.ref  NONE;

  1406 	val len= Unsynchronized.ref  0;

  1407 	val vl= Unsynchronized.ref  [];

  1408 	val vh= Unsynchronized.ref  [];

  1409 	val vtemp= Unsynchronized.ref  [];

  1410 	val i= Unsynchronized.ref  0;

  1411     in

  1412     (

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

  1414 	 (

  1415 	  len:=length(v);

  1416 	  vh:=v;

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

  1419 	      (

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

  1421 		   (

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

  1423 		    )

  1424 	       else

  1425 		   (

  1426 		    vl:=0::(!vl)

  1427 		    );

  1428 		   vh:=tl(!vh);

  1429 		   i:=(!i)+1

  1430 		   );

  1431 	      SOME [(1,rev(!vl))] handle _ => NONE

  1432 	      )

  1433      else error ("RATIONALS_TERM2COEF_EXCEPTION 1: Invalid term")

  1434 	 )

  1435     end

  1436   | term2coef' (Const ("Groups.plus_class.plus",_) $ t1 $ t2) v :mv_poly option=

  1437     (

  1438      SOME ((the(term2coef' t1 v)) @ (the(term2coef' t2 v))) handle _ => NONE

  1439 	 )

  1440   | term2coef' (Const ("Groups.minus_class.minus",_) $ t1 $ t2) v :mv_poly option=

  1441     (

  1442      SOME ((the(term2coef' t1 v)) @ mv_skalar_mul((the(term2coef' t2 v)),1)) handle _ => NONE

  1443 	 )

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

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

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

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

  1449     else false;

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

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

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

  1455 fun mk_monom v' p vs =

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

  1457     in map (conv p) vs end;

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

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

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

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

  1463     if is_numeral str

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

  1465     else SOME [(~1, mk_monom str 1 v)]

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

  1468     let

  1469 	val x= Unsynchronized.ref  NONE;

  1470 	val len= Unsynchronized.ref  0;

  1471 	val vl= Unsynchronized.ref  [];

  1472 	val vh= Unsynchronized.ref  [];

  1473 	val i= Unsynchronized.ref  0;

  1474     in

  1475 	if is_numeral str then

  1476 	    (

  1477 	     SOME [(((the o int_of_str) str),mv_null2 v)] handle _ => NONE

  1478 		 )

  1479 	else (* variable *)

  1480 	    (

  1481 	     len:=length v;

  1482 	     vh:= v;

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

  1484 		 (

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

  1486 		      (

  1487 		       vl:=1::(!vl)

  1488 		       )

  1489 		  else

  1490 		      (

  1491 		       vl:=0::(!vl)

  1492 		       );

  1493 		      vh:=tl(!vh);

  1494 		      i:=(!i)+1

  1495 		      );

  1496 		 SOME [(1,rev(!vl))] handle _ => NONE

  1497 	    )

  1498     end

  1499   | term2poly' (Const ("Groups.times_class.times",_) $ t1 $ t2) v :mv_poly option=

  1500     let

  1501 	val t1pp= Unsynchronized.ref  [];

  1502 	val t2pp= Unsynchronized.ref  [];

  1503 	val t1c= Unsynchronized.ref  0;

  1504 	val t2c= Unsynchronized.ref  0;

  1505     in

  1506 	(

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

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

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

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

  1512 	 SOME [( (!t1c)*(!t2c) ,( (map op+) ((!t1pp)~~(!t2pp)) ) )]

  1513 	 handle _ => NONE

  1515 	 )

  1516     end

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

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

  1519     let

  1520 	val x= Unsynchronized.ref  NONE;

  1521 	val len= Unsynchronized.ref  0;

  1522 	val vl= Unsynchronized.ref  [];

  1523 	val vh= Unsynchronized.ref  [];

  1524 	val vtemp= Unsynchronized.ref  [];

  1525 	val i= Unsynchronized.ref  0;

  1526     in

  1527     (

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

  1529 	 (

  1530 	  len:=length(v);

  1531 	  vh:=v;

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

  1534 	      (

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

  1536 		   (

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

  1538 		    )

  1539 	       else

  1540 		   (

  1541 		    vl:=0::(!vl)

  1542 		    );

  1543 		   vh:=tl(!vh);

  1544 		   i:=(!i)+1

  1545 		   );

  1546 	      SOME [(1,rev(!vl))] handle _ => NONE

  1547 	      )

  1548      else error ("RATIONALS_TERM2POLY_EXCEPTION 1: Invalid term")

  1549 	 )

  1550     end

  1551   | term2poly' (Const ("Groups.plus_class.plus",_) $ t1 $ t2) v :mv_poly option =

  1552     (

  1553      SOME ((the(term2poly' t1 v)) @ (the(term2poly' t2 v))) handle _ => NONE

  1554 	 )

  1555   | term2poly' (Const ("Groups.minus_class.minus",_) $ t1 $ t2) v :mv_poly option =

  1556     (

  1557      SOME ((the(term2poly' t1 v)) @ mv_skalar_mul((the(term2poly' t2 v)),~1)) handle _ => NONE

  1558 	 )

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

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

  1562     term2poly

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

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

  1565     			       DIE DU MIR LETZTES MAL GEGEBEN HAST*)

  1566     	 mv_monom list        (*internal representation           *)

  1567     		  option      (*the translation may fail with NONE*)

  1568 .*)

  1569 fun term2poly (t:term) v =

  1570      if is_polynomial t then term2poly' t v

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

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

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

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

  1577 fun expanded2poly (t:term) v =

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

  1579     (*else error ("RATIONALS_EXPANDED2POLY_EXCEPTION: Invalid Polynomial")*);

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

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

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

  1585 fun powerproduct2term(xs,v) =

  1586     let

  1587 	val xss= Unsynchronized.ref  xs;

  1588 	val vv= Unsynchronized.ref  v;

  1589     in

  1590 	(

  1591 	 while hd(!xss)=0 do

  1592 	     (

  1593 	      xss:=tl(!xss);

  1594 	      vv:=tl(!vv)

  1595 	      );

  1597 	 if list_is_null(tl(!xss)) then

  1598 	     (

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

  1600 	      else

  1601 		  (

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

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

  1604 		   )

  1605 	      )

  1606 	 else

  1607 	     (

  1608 	      if hd(!xss)=1 then

  1609 		  (

  1610 		   Const("Groups.times_class.times",[HOLogic.realT,HOLogic.realT]--->HOLogic.realT) $

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

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

  1613 		   )

  1614 	      else

  1615 		  (

  1616 		   Const("Groups.times_class.times",[HOLogic.realT,HOLogic.realT]--->HOLogic.realT) $

  1617 		   (

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

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

  1620 		    ) $

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

  1622 		   )

  1623 	      )

  1624 	 )

  1625     end;

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

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

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

  1630     else

  1631 	(

  1632 	 if list_is_null(e) then

  1633 	     (

  1634 	      Free(str_of_int c,HOLogic.realT)

  1635 	      )

  1636 	 else

  1637 	     (

  1638 	      if c=1 then

  1639 		  (

  1640 		   powerproduct2term(e,v)

  1641 		   )

  1642 	      else

  1643 		  (

  1644 		   Const("Groups.times_class.times",[HOLogic.realT,HOLogic.realT]--->HOLogic.realT) $

  1645 		   Free(str_of_int c,HOLogic.realT)  $

  1646 		   powerproduct2term(e,v)

  1647 		   )

  1648 		  )

  1649 	     );*)

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

  1653     if list_is_null is

  1654     then

  1655 	if i >= 0

  1656 	then Free (str_of_int i, HOLogic.realT)

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

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

  1659     else

  1660 	if i > 0

  1661 	then Const ("Groups.times_class.times", [HOLogic.realT,HOLogic.realT]---> HOLogic.realT) $

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

  1663 		   powerproduct2term(is, v)

  1664 	else Const ("Groups.times_class.times", [HOLogic.realT,HOLogic.realT]---> HOLogic.realT) $

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

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

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

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

  1669     if list_is_null is

  1670     then Free (str_of_int i, HOLogic.realT)

  1671     else if i = 1

  1672     then powerproduct2term (is, vs)

  1673     else Const ("Groups.times_class.times", [HOLogic.realT, HOLogic.realT] ---> HOLogic.realT) $

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

  1675 	       powerproduct2term (is, vs);

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

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

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

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

  1681     Const("Groups.plus_class.plus", [HOLogic.realT, HOLogic.realT] ---> HOLogic.realT) $

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

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

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

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

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

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

  1690     else

  1691 	(

  1692 	 if list_is_null(e) then

  1693 	     (

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

  1695 	      )

  1696 	 else

  1697 	     (

  1698 	      if abs(c)=1 then

  1699 		  (

  1700 		   powerproduct2term(e,v)

  1701 		   )

  1702 	      else

  1703 		  (

  1704 		   Const("Groups.times_class.times",[HOLogic.realT,HOLogic.realT]--->HOLogic.realT) $

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

  1706 		   powerproduct2term(e,v)

  1707 		   )

  1708 		  )

  1709 	     );

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

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

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

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

  1715     if c1<0 then

  1716 	Const("Groups.minus_class.minus",[HOLogic.realT,HOLogic.realT]--->HOLogic.realT) $

  1717 	exp2term'(others,vars) $

  1718 	(

  1719 	 monom2term2((c1,e1),vars)

  1720 	 )

  1721     else

  1722 	Const("Groups.plus_class.plus",[HOLogic.realT,HOLogic.realT]--->HOLogic.realT) $

  1723 	exp2term'(others,vars) $

  1724 	(

  1725 	 monom2term2((c1,e1),vars)

  1726 	 );

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

  1729     a term in polynomial representation .*)

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

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

  1733 fun polynomial2expanded t =

  1734     (let

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

  1736     in

  1737 	SOME (poly2expanded (the (term2poly t vars), vars))

  1738     end) handle _ => NONE;

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

  1741 fun expanded2polynomial t =

  1742     (let

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

  1744     in

  1745 	SOME (poly2term (the (expanded2poly t vars), vars))

  1746     end) handle _ => NONE;

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

  1750 fun step_cancel (t as Const ("Rings.inverse_class.divide",_) $ p1 $ p2) =

  1751     let

  1752 	val p1' = Unsynchronized.ref [];

  1753 	val p2' = Unsynchronized.ref [];

  1754 	val p3  = Unsynchronized.ref []

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

  1756     in

  1757 	(

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

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

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

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

  1762 	     (

  1763 	      Const ("Rings.inverse_class.divide",[HOLogic.realT,HOLogic.realT]--->HOLogic.realT) $ p1 $ p2

  1764 	      )

  1765 	 else

  1766 	     (

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

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

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

  1772 	      (

  1773 	       Const ("Rings.inverse_class.divide",[HOLogic.realT,HOLogic.realT]--->HOLogic.realT)

  1774 	       $

  1775 	       (

  1776 		Const ("Groups.times_class.times",[HOLogic.realT,HOLogic.realT]--->HOLogic.realT) $

  1777 		poly2term(!p1',vars) $

  1778 		poly2term(!p3,vars)

  1779 		)

  1780 	       $

  1781 	       (

  1782 		Const ("Groups.times_class.times",[HOLogic.realT,HOLogic.realT]--->HOLogic.realT) $

  1783 		poly2term(!p2',vars) $

  1784 		poly2term(!p3,vars)

  1785 		)

  1786 	       )

  1787 	      else

  1788 	      (

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

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

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

  1792 	       (

  1793 		Const ("Rings.inverse_class.divide",[HOLogic.realT,HOLogic.realT]--->HOLogic.realT)

  1794 		$

  1795 		(

  1796 		 Const ("Groups.times_class.times",[HOLogic.realT,HOLogic.realT]--->HOLogic.realT) $

  1797 		 poly2term(!p1',vars) $

  1798 		 poly2term(!p3,vars)

  1799 		 )

  1800 		$

  1801 		(

  1802 		 Const ("Groups.times_class.times",[HOLogic.realT,HOLogic.realT]--->HOLogic.realT) $

  1803 		 poly2term(!p2',vars) $

  1804 		 poly2term(!p3,vars)

  1805 		 )

  1806 		)

  1807 	       )

  1808 	      )

  1809 	     )

  1810     end

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

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

  1815 fun step_cancel_expanded (t as Const ("Rings.inverse_class.divide",_) $ p1 $ p2) =

  1816     let

  1817 	val p1' = Unsynchronized.ref [];

  1818 	val p2' = Unsynchronized.ref [];

  1819 	val p3  = Unsynchronized.ref []

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

  1821     in

  1822 	(

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

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

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

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

  1827 	     (

  1828 	      Const ("Rings.inverse_class.divide",[HOLogic.realT,HOLogic.realT]--->HOLogic.realT) $ p1 $ p2

  1829 	      )

  1830 	 else

  1831 	     (

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

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

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

  1837 	      (

  1838 	       Const ("Rings.inverse_class.divide",[HOLogic.realT,HOLogic.realT]--->HOLogic.realT)

  1839 	       $

  1840 	       (

  1841 		Const ("Groups.times_class.times",[HOLogic.realT,HOLogic.realT]--->HOLogic.realT) $

  1842 		poly2expanded(!p1',vars) $

  1843 		poly2expanded(!p3,vars)

  1844 		)

  1845 	       $

  1846 	       (

  1847 		Const ("Groups.times_class.times",[HOLogic.realT,HOLogic.realT]--->HOLogic.realT) $

  1848 		poly2expanded(!p2',vars) $

  1849 		poly2expanded(!p3,vars)

  1850 		)

  1851 	       )

  1852 	      else

  1853 	      (

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

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

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

  1857 	       (

  1858 		Const ("Rings.inverse_class.divide",[HOLogic.realT,HOLogic.realT]--->HOLogic.realT)

  1859 		$

  1860 		(

  1861 		 Const ("Groups.times_class.times",[HOLogic.realT,HOLogic.realT]--->HOLogic.realT) $

  1862 		 poly2expanded(!p1',vars) $

  1863 		 poly2expanded(!p3,vars)

  1864 		 )

  1865 		$

  1866 		(

  1867 		 Const ("Groups.times_class.times",[HOLogic.realT,HOLogic.realT]--->HOLogic.realT) $

  1868 		 poly2expanded(!p2',vars) $

  1869 		 poly2expanded(!p3,vars)

  1870 		 )

  1871 		)

  1872 	       )

  1873 	      )

  1874 	     )

  1875     end

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

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

  1879 fun direct_cancel (t as Const ("Rings.inverse_class.divide",_) $ p1 $ p2) =

  1880     let

  1881 	val p1' = Unsynchronized.ref [];

  1882 	val p2' = Unsynchronized.ref [];

  1883 	val p3  = Unsynchronized.ref []

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

  1885     in

  1886 	(

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

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

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

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

  1892 	     (

  1893 	      (Const ("Rings.inverse_class.divide",[HOLogic.realT,HOLogic.realT]--->HOLogic.realT) $ p1 $ p2,[])

  1894 	      )

  1895 	 else

  1896 	     (

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

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

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

  1900 	      (

  1901 	       (

  1902 		Const ("Rings.inverse_class.divide",[HOLogic.realT,HOLogic.realT]--->HOLogic.realT)

  1903 		$

  1904 		(

  1905 		 poly2term((!p1'),vars)

  1906 		 )

  1907 		$

  1908 		(

  1909 		 poly2term((!p2'),vars)

  1910 		 )

  1911 		)

  1912 	       ,

  1913 	       if mv_grad(!p3)>0 then

  1914 		   [

  1915 		    (

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

  1917 		     (

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

  1919 		      poly2term((!p3),vars) $

  1920 		      Free("0",HOLogic.realT)

  1921 		      )

  1922 		     )

  1923 		    ]

  1924 	       else

  1925 		   []

  1926 		   )

  1927 	      else

  1928 		  (

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

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

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

  1932 		       (

  1933 			Const ("Rings.inverse_class.divide",[HOLogic.realT,HOLogic.realT]--->HOLogic.realT)

  1934 			$

  1935 			(

  1936 			 poly2term((!p1'),vars)

  1937 			 )

  1938 			$

  1939 			(

  1940 			 poly2term((!p2'),vars)

  1941 			 )

  1942 			,

  1943 			if mv_grad(!p3)>0 then

  1944 			    [

  1945 			     (

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

  1947 			      (

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

  1949 			       poly2term((!p3),vars) $

  1950 			       Free("0",HOLogic.realT)

  1951 			       )

  1952 			      )

  1953 			     ]

  1954 			else

  1955 			    []

  1956 			    )

  1957 		       )

  1958 		  )

  1959 	     )

  1960     end

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

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

  1964 fun direct_cancel_expanded (t as Const ("Rings.inverse_class.divide",_) $ p1 $ p2) =

  1965     let

  1966 	val p1' = Unsynchronized.ref [];

  1967 	val p2' = Unsynchronized.ref [];

  1968 	val p3  = Unsynchronized.ref []

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

  1970     in

  1971 	(

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

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

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

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

  1977 	     (

  1978 	      (Const ("Rings.inverse_class.divide",[HOLogic.realT,HOLogic.realT]--->HOLogic.realT) $ p1 $ p2,[])

  1979 	      )

  1980 	 else

  1981 	     (

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

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

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

  1985 	      (

  1986 	       (

  1987 		Const ("Rings.inverse_class.divide",[HOLogic.realT,HOLogic.realT]--->HOLogic.realT)

  1988 		$

  1989 		(

  1990 		 poly2expanded((!p1'),vars)

  1991 		 )

  1992 		$

  1993 		(

  1994 		 poly2expanded((!p2'),vars)

  1995 		 )

  1996 		)

  1997 	       ,

  1998 	       if mv_grad(!p3)>0 then

  1999 		   [

  2000 		    (

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

  2002 		     (

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

  2004 		      poly2expanded((!p3),vars) $

  2005 		      Free("0",HOLogic.realT)

  2006 		      )

  2007 		     )

  2008 		    ]

  2009 	       else

  2010 		   []

  2011 		   )

  2012 	      else

  2013 		  (

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

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

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

  2017 		       (

  2018 			Const ("Rings.inverse_class.divide",[HOLogic.realT,HOLogic.realT]--->HOLogic.realT)

  2019 			$

  2020 			(

  2021 			 poly2expanded((!p1'),vars)

  2022 			 )

  2023 			$

  2024 			(

  2025 			 poly2expanded((!p2'),vars)

  2026 			 )

  2027 			,

  2028 			if mv_grad(!p3)>0 then

  2029 			    [

  2030 			     (

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

  2032 			      (

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

  2034 			       poly2expanded((!p3),vars) $

  2035 			       Free("0",HOLogic.realT)

  2036 			       )

  2037 			      )

  2038 			     ]

  2039 			else

  2040 			    []

  2041 			    )

  2042 		       )

  2043 		  )

  2044 	     )

  2045     end

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

  2049 (*. adds two fractions .*)

  2050 fun add_fract ((Const("Rings.inverse_class.divide",_) $ t11 $ t12),(Const("Rings.inverse_class.divide",_) $ t21 $ t22)) =

  2051     let

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

  2053 	val t11'= Unsynchronized.ref  (the(term2poly t11 vars));

  2054 val _= tracing"### add_fract: done t11"

  2055 	val t12'= Unsynchronized.ref  (the(term2poly t12 vars));

  2056 val _= tracing"### add_fract: done t12"

  2057 	val t21'= Unsynchronized.ref  (the(term2poly t21 vars));

  2058 val _= tracing"### add_fract: done t21"

  2059 	val t22'= Unsynchronized.ref  (the(term2poly t22 vars));

  2060 val _= tracing"### add_fract: done t22"

  2061 	val den= Unsynchronized.ref  [];

  2062 	val nom= Unsynchronized.ref  [];

  2063 	val m1= Unsynchronized.ref  [];

  2064 	val m2= Unsynchronized.ref  [];

  2065     in

  2067 	(

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

  2069 tracing"### add_fract: done sort mv_lcm";

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

  2071 tracing"### add_fract: done sort mv_division t12";

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

  2073 tracing"### add_fract: done sort mv_division t22";

  2074 	 nom :=sort (mv_geq LEX_)

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

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

  2077 				       LEX_),

  2078 				LEX_));

  2079 tracing"### add_fract: done sort mv_add";

  2080 	 (

  2081 	  Const ("Rings.inverse_class.divide",[HOLogic.realT,HOLogic.realT]--->HOLogic.realT)

  2082 	  $

  2083 	  (

  2084 	   poly2term((!nom),vars)

  2085 	   )

  2086 	  $

  2087 	  (

  2088 	   poly2term((!den),vars)

  2089 	   )

  2090 	  )

  2091 	 )

  2092     end

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

  2095 (*. adds two expanded fractions .*)

  2096 fun add_fract_exp ((Const("Rings.inverse_class.divide",_) $ t11 $ t12),(Const("Rings.inverse_class.divide",_) $ t21 $ t22)) =

  2097     let

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

  2099 	val t11'= Unsynchronized.ref  (the(expanded2poly t11 vars));

  2100 	val t12'= Unsynchronized.ref  (the(expanded2poly t12 vars));

  2101 	val t21'= Unsynchronized.ref  (the(expanded2poly t21 vars));

  2102 	val t22'= Unsynchronized.ref  (the(expanded2poly t22 vars));

  2103 	val den= Unsynchronized.ref  [];

  2104 	val nom= Unsynchronized.ref  [];

  2105 	val m1= Unsynchronized.ref  [];

  2106 	val m2= Unsynchronized.ref  [];

  2107     in

  2109 	(

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

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

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

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

  2114 	 (

  2115 	  Const ("Rings.inverse_class.divide",[HOLogic.realT,HOLogic.realT]--->HOLogic.realT)

  2116 	  $

  2117 	  (

  2118 	   poly2expanded((!nom),vars)

  2119 	   )

  2120 	  $

  2121 	  (

  2122 	   poly2expanded((!den),vars)

  2123 	   )

  2124 	  )

  2125 	 )

  2126     end

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

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

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

  2131   | add_list_of_fractions [x]= direct_cancel x

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

  2133     let

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

  2135 val _= tracing"### add_list_of_fractions xs: has done direct_cancel(x)";

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

  2137 val _= tracing"### add_list_of_fractions xs: has done direct_cancel(y)";

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

  2139 val _= tracing"### add_list_of_fractions xs: has done add_list_of_fraction xs";

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

  2141 val _= tracing"### add_list_of_fractions xs: has done direct_cancel(t3a)";

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

  2143     in

  2144 	(tracing"### add_list_of_fractions in";

  2145 	 (

  2146 	 (t4a,rest)

  2147 	 )

  2148 	 )

  2149     end;

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

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

  2153   | add_list_of_fractions_exp [x]= direct_cancel_expanded x

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

  2155     let

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

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

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

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

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

  2161     in

  2162 	(

  2163 	 (t4a,rest)

  2164 	 )

  2165     end;

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

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

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

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

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

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

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

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

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

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

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

  2181   | termlist2denominators xs =

  2182     let

  2183 	val xxs= Unsynchronized.ref  xs;

  2184 	val var= Unsynchronized.ref  [];

  2185     in

  2186 	var:=[];

  2187 	while length(!xxs)>0 do

  2188 	    (

  2189 	     let

  2190 		 val (t as Const ("Rings.inverse_class.divide",_) $ p1x $ p2x)=hd(!xxs);

  2191 	     in

  2192 		 (

  2193 		  xxs:=tl(!xxs);

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

  2195 		  )

  2196 	     end

  2197 	     );

  2198 	    (t2d(xs,!var),!var)

  2199     end;

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

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

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

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

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

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

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

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

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

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

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

  2215   | termlist2denominators xs =

  2216     let

  2217 	val xxs= Unsynchronized.ref  xs;

  2218 	val var= Unsynchronized.ref  [];

  2219     in

  2220 	var:=[];

  2221 	while length(!xxs)>0 do

  2222 	    (

  2223 	     let

  2224 		 val (t as Const ("Rings.inverse_class.divide",_) $ p1x $ p2x)=hd(!xxs);

  2225 	     in

  2226 		 (

  2227 		  xxs:=tl(!xxs);

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

  2229 		  )

  2230 	     end

  2231 	     );

  2232 	    (t2d(xs,!var),!var)

  2233     end;

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

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

  2237   | termlist2denominators_exp xs =

  2238     let

  2239 	val xxs= Unsynchronized.ref  xs;

  2240 	val var= Unsynchronized.ref  [];

  2241     in

  2242 	var:=[];

  2243 	while length(!xxs)>0 do

  2244 	    (

  2245 	     let

  2246 		 val (t as Const ("Rings.inverse_class.divide",_) $ p1x $ p2x)=hd(!xxs);

  2247 	     in

  2248 		 (

  2249 		  xxs:=tl(!xxs);

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

  2251 		  )

  2252 	     end

  2253 	     );

  2254 	    (t2d_exp(xs,!var),!var)

  2255     end;

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

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

  2259     let

  2260 	val (t as Const ("Rings.inverse_class.divide",_) $ p1' $ p2')=x;

  2261 	val p2= sort (mv_geq LEX_) (the(term2poly 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 ("Groups.plus_class.plus",[HOLogic.realT,HOLogic.realT]--->HOLogic.realT)

  2269 		 $

  2270 		 (

  2271 		  Const ("Rings.inverse_class.divide",[HOLogic.realT,HOLogic.realT]--->HOLogic.realT)

  2272 		  $

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

  2274 		  $

  2275 		  den

  2276 		  )

  2277 		 $

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

  2279 		,

  2280 		[]

  2281 		)

  2282 	    else

  2283 		(

  2284 		 Const ("Groups.plus_class.plus",[HOLogic.realT,HOLogic.realT]--->HOLogic.realT)

  2285 		 $

  2286 		 (

  2287 		  Const ("Rings.inverse_class.divide",[HOLogic.realT,HOLogic.realT]--->HOLogic.realT)

  2288 		  $

  2289 		  (

  2290 		   Const ("Groups.times_class.times",[HOLogic.realT,HOLogic.realT]--->HOLogic.realT) $

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

  2292 		   poly2term(p3,var)

  2293 		   )

  2294 		  $

  2295 		  (

  2296 		   den

  2297 		   )

  2298 		  )

  2299 		 $

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

  2301 		,

  2302 		[]

  2303 		)

  2304 	else

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

  2306 		(

  2307 		 (

  2308 		  Const ("Rings.inverse_class.divide",[HOLogic.realT,HOLogic.realT]--->HOLogic.realT)

  2309 		  $

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

  2311 		  $

  2312 		  den

  2313 		  )

  2314 		 ,

  2315 		 []

  2316 		 )

  2317 	     else

  2318 		 (

  2319 		  Const ("Rings.inverse_class.divide",[HOLogic.realT,HOLogic.realT]--->HOLogic.realT)

  2320 		  $

  2321 		  (

  2322 		   Const ("Groups.times_class.times",[HOLogic.realT,HOLogic.realT]--->HOLogic.realT) $

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

  2324 		   poly2term(p3,var)

  2325 		   )

  2326 		  $

  2327 		  den

  2328 		  ,

  2329 		  []

  2330 		  )

  2331     end;

  2333 (*. same as com_den, this time for expanded forms .*)

  2334 fun com_den_exp(x::xs,denom,den,var)=

  2335     let

  2336 	val (t as Const ("Rings.inverse_class.divide",_) $ p1' $ p2')=x;

  2337 	val p2= sort (mv_geq LEX_) (the(expanded2poly p2' var));

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

  2339 	val p1var=get_vars(p1');

  2340     in

  2341 	if length(xs)>0 then

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

  2343 		(

  2344 		 Const ("Groups.plus_class.plus",[HOLogic.realT,HOLogic.realT]--->HOLogic.realT)

  2345 		 $

  2346 		 (

  2347 		  Const ("Rings.inverse_class.divide",[HOLogic.realT,HOLogic.realT]--->HOLogic.realT)

  2348 		  $

  2349 		  poly2expanded(the(expanded2poly p1' p1var),p1var)

  2350 		  $

  2351 		  den

  2352 		  )

  2353 		 $

  2354 		 #1(com_den_exp(xs,denom,den,var))

  2355 		,

  2356 		[]

  2357 		)

  2358 	    else

  2359 		(

  2360 		 Const ("Groups.plus_class.plus",[HOLogic.realT,HOLogic.realT]--->HOLogic.realT)

  2361 		 $

  2362 		 (

  2363 		  Const ("Rings.inverse_class.divide",[HOLogic.realT,HOLogic.realT]--->HOLogic.realT)

  2364 		  $

  2365 		  (

  2366 		   Const ("Groups.times_class.times",[HOLogic.realT,HOLogic.realT]--->HOLogic.realT) $

  2367 		   poly2expanded(the(expanded2poly p1' p1var),p1var) $

  2368 		   poly2expanded(p3,var)

  2369 		   )

  2370 		  $

  2371 		  (

  2372 		   den

  2373 		   )

  2374 		  )

  2375 		 $

  2376 		 #1(com_den_exp(xs,denom,den,var))

  2377 		,

  2378 		[]

  2379 		)

  2380 	else

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

  2382 		(

  2383 		 (

  2384 		  Const ("Rings.inverse_class.divide",[HOLogic.realT,HOLogic.realT]--->HOLogic.realT)

  2385 		  $

  2386 		  poly2expanded(the(expanded2poly p1' p1var),p1var)

  2387 		  $

  2388 		  den

  2389 		  )

  2390 		 ,

  2391 		 []

  2392 		 )

  2393 	     else

  2394 		 (

  2395 		  Const ("Rings.inverse_class.divide",[HOLogic.realT,HOLogic.realT]--->HOLogic.realT)

  2396 		  $

  2397 		  (

  2398 		   Const ("Groups.times_class.times",[HOLogic.realT,HOLogic.realT]--->HOLogic.realT) $

  2399 		   poly2expanded(the(expanded2poly p1' p1var),p1var) $

  2400 		   poly2expanded(p3,var)

  2401 		   )

  2402 		  $

  2403 		  den

  2404 		  ,

  2405 		  []

  2406 		  )

  2407     end;

  2409 (* wird aktuell nicht mehr gebraucht, bei rückänderung schon

  2410 -------------------------------------------------------------

  2411 (* WN0210???SK brauch ma des überhaupt *)

  2412 fun com_den2(x::xs,denom,den,var)=

  2413     let

  2414 	val (t as Const ("Rings.inverse_class.divide",_) $ p1' $ p2')=x;

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

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

  2417 	val p1var=get_vars(p1');

  2418     in

  2419 	if length(xs)>0 then

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

  2421 		(

  2422 		 Const ("Groups.plus_class.plus",[HOLogic.realT,HOLogic.realT]--->HOLogic.realT) $

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

  2424 		 com_den2(xs,denom,den,var)

  2425 		)

  2426 	    else

  2427 		(

  2428 		 Const ("Groups.plus_class.plus",[HOLogic.realT,HOLogic.realT]--->HOLogic.realT) $

  2429 		 (

  2430 		   let

  2431 		       val p3'=poly2term(p3,var);

  2432 		       val vars= (((map free2str) o vars) p1') union (((map free2str) o vars) p3');

  2433 		   in

  2434 		       poly2term(sort (mv_geq LEX_) (mv_mul(the(term2poly p1' vars) ,the(term2poly p3' vars),LEX_)),vars)

  2435 		   end

  2436 		  ) $

  2437 		 com_den2(xs,denom,den,var)

  2438 		)

  2439 	else

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

  2441 		(

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

  2443 		 )

  2444 	     else

  2445 		 (

  2446 		   let

  2447 		       val p3'=poly2term(p3,var);

  2448 		       val vars= (((map free2str) o vars) p1') union (((map free2str) o vars) p3');

  2449 		   in

  2450 		       poly2term(sort (mv_geq LEX_) (mv_mul(the(term2poly p1' vars) ,the(term2poly p3' vars),LEX_)),vars)

  2451 		   end

  2452 		  )

  2453     end;

  2455 (* WN0210???SK brauch ma des überhaupt *)

  2456 fun com_den_exp2(x::xs,denom,den,var)=

  2457     let

  2458 	val (t as Const ("Rings.inverse_class.divide",_) $ p1' $ p2')=x;

  2459 	val p2= sort (mv_geq LEX_) (the(expanded2poly p2' var));

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

  2461 	val p1var=get_vars p1';

  2462     in

  2463 	if length(xs)>0 then

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

  2465 		(

  2466 		 Const ("Groups.plus_class.plus",[HOLogic.realT,HOLogic.realT]--->HOLogic.realT) $

  2467 		 poly2expanded(the (expanded2poly p1' p1var),p1var) $

  2468 		 com_den_exp2(xs,denom,den,var)

  2469 		)

  2470 	    else

  2471 		(

  2472 		 Const ("Groups.plus_class.plus",[HOLogic.realT,HOLogic.realT]--->HOLogic.realT) $

  2473 		 (

  2474 		   let

  2475 		       val p3'=poly2expanded(p3,var);

  2476 		       val vars= (((map free2str) o vars) p1') union (((map free2str) o vars) p3');

  2477 		   in

  2478 		       poly2expanded(sort (mv_geq LEX_) (mv_mul(the(expanded2poly p1' vars) ,the(expanded2poly p3' vars),LEX_)),vars)

  2479 		   end

  2480 		  ) $

  2481 		 com_den_exp2(xs,denom,den,var)

  2482 		)

  2483 	else

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

  2485 		(

  2486 		 poly2expanded(the (expanded2poly p1' p1var),p1var)

  2487 		 )

  2488 	     else

  2489 		 (

  2490 		   let

  2491 		       val p3'=poly2expanded(p3,var);

  2492 		       val vars= (((map free2str) o vars) p1') union (((map free2str) o vars) p3');

  2493 		   in

  2494 		       poly2expanded(sort (mv_geq LEX_) (mv_mul(the(expanded2poly p1' vars) ,the(expanded2poly p3' vars),LEX_)),vars)

  2495 		   end

  2496 		  )

  2497     end;

  2498 ---------------------------------------------------------*)

  2501 (*. searches for an element y of a list ys, which has an gcd not 1 with x .*)

  2502 fun exists_gcd (x,[]) = false

  2503   | exists_gcd (x,y::ys) = if mv_gcd x y = [(1,mv_null2(#2(hd(x))))] then  exists_gcd (x,ys)

  2504 			   else true;

  2506 (*. divides each element of the list xs with y .*)

  2507 fun list_div ([],y) = []

  2508   | list_div (x::xs,y) =

  2509     let

  2510 	val (d,r)=mv_division(x,y,LEX_);

  2511     in

  2512 	if r=[] then

  2513 	    d::list_div(xs,y)

  2514 	else x::list_div(xs,y)

  2515     end;

  2517 (*. checks if x is in the list ys .*)

  2518 fun in_list (x,[]) = false

  2519   | in_list (x,y::ys) = if x=y then true

  2520 			else in_list(x,ys);

  2522 (*. deletes all equal elements of the list xs .*)

  2523 fun kill_equal [] = []

  2524   | kill_equal (x::xs) = if in_list(x,xs) orelse x=[(1,mv_null2(#2(hd(x))))] then kill_equal(xs)

  2525 			 else x::kill_equal(xs);

  2527 (*. searches for new factors .*)

  2528 fun new_factors [] = []

  2529   | new_factors (list:mv_poly list):mv_poly list =

  2530     let

  2531 	val l = kill_equal list;

  2532 	val len = length(l);

  2533     in

  2534 	if len>=2 then

  2535 	    (

  2536 	     let

  2537 		 val x::y::xs=l;

  2538 		 val gcd=mv_gcd x y;

  2539 	     in

  2540 		 if gcd=[(1,mv_null2(#2(hd(x))))] then

  2541 		     (

  2542 		      if exists_gcd(x,xs) then new_factors (y::xs @ [x])

  2543 		      else x::new_factors(y::xs)

  2544 	             )

  2545 		 else gcd::new_factors(kill_equal(list_div(x::y::xs,gcd)))

  2546 	     end

  2547 	     )

  2548 	else

  2549 	    if len=1 then [hd(l)]

  2550 	    else []

  2551     end;

  2553 (*. gets the factors of a list .*)

  2554 fun get_factors x = new_factors x;

  2556 (*. multiplies the elements of the list .*)

  2557 fun multi_list [] = []

  2558   | multi_list (x::xs) = if xs=[] then x

  2559 			 else mv_mul(x,multi_list xs,LEX_);

  2561 (*. makes a term out of the elements of the list (polynomial representation) .*)

  2562 fun make_term ([],vars) = Free(str_of_int 0,HOLogic.realT)

  2563   | make_term ((x::xs),vars) = if length(xs)=0 then poly2term(sort (mv_geq LEX_) (x),vars)

  2564 			       else

  2565 				   (

  2566 				    Const ("Groups.times_class.times",[HOLogic.realT,HOLogic.realT]--->HOLogic.realT) $

  2567 				    poly2term(sort (mv_geq LEX_) (x),vars) $

  2568 				    make_term(xs,vars)

  2569 				    );

  2571 (*. factorizes the denominator (polynomial representation) .*)

  2572 fun factorize_den (l,den,vars) =

  2573     let

  2574 	val factor_list=kill_equal( (get_factors l));

  2575 	val mlist=multi_list(factor_list);

  2576 	val (last,rest)=mv_division(den,multi_list(factor_list),LEX_);

  2577     in

  2578 	if rest=[] then

  2579 	    (

  2580 	     if last=[(1,mv_null2(vars))] then make_term(factor_list,vars)

  2581 	     else make_term(last::factor_list,vars)

  2582 	     )

  2583 	else error ("RATIONALS_FACTORIZE_DEN_EXCEPTION: Invalid factor by division")

  2584     end;

  2586 (*. makes a term out of the elements of the list (expanded polynomial representation) .*)

  2587 fun make_exp ([],vars) = Free(str_of_int 0,HOLogic.realT)

  2588   | make_exp ((x::xs),vars) = if length(xs)=0 then poly2expanded(sort (mv_geq LEX_) (x),vars)

  2589 			       else

  2590 				   (

  2591 				    Const ("Groups.times_class.times",[HOLogic.realT,HOLogic.realT]--->HOLogic.realT) $

  2592 				    poly2expanded(sort (mv_geq LEX_) (x),vars) $

  2593 				    make_exp(xs,vars)

  2594 				    );

  2596 (*. factorizes the denominator (expanded polynomial representation) .*)

  2597 fun factorize_den_exp (l,den,vars) =

  2598     let

  2599 	val factor_list=kill_equal( (get_factors l));

  2600 	val mlist=multi_list(factor_list);

  2601 	val (last,rest)=mv_division(den,multi_list(factor_list),LEX_);

  2602     in

  2603 	if rest=[] then

  2604 	    (

  2605 	     if last=[(1,mv_null2(vars))] then make_exp(factor_list,vars)

  2606 	     else make_exp(last::factor_list,vars)

  2607 	     )

  2608 	else error ("RATIONALS_FACTORIZE_DEN_EXP_EXCEPTION: Invalid factor by division")

  2609     end;

  2611 (*. calculates the common denominator of all elements of the list and multiplies .*)

  2612 (*. the nominators and denominators with the correct factor .*)

  2613 (*. (polynomial representation) .*)

  2614 fun step_add_list_of_fractions []=(Free("0",HOLogic.realT),[]:term list)

  2615   | step_add_list_of_fractions [x]= error ("RATIONALS_STEP_ADD_LIST_OF_FRACTIONS_EXCEPTION: Nothing to add")

  2616   | step_add_list_of_fractions (xs) =

  2617     let

  2618         val den_list=termlist2denominators (xs); (* list of denominators *)

  2619 	val (denom,var)=calc_lcm(den_list);      (* common denominator *)

  2620 	val den=factorize_den(#1(den_list),denom,var); (* faktorisierter Nenner !!! *)

  2621     in

  2622 	com_den(xs,denom,den,var)

  2623     end;

  2625 (*. calculates the common denominator of all elements of the list and multiplies .*)

  2626 (*. the nominators and denominators with the correct factor .*)

  2627 (*. (expanded polynomial representation) .*)

  2628 fun step_add_list_of_fractions_exp []  = (Free("0",HOLogic.realT),[]:term list)

  2629   | step_add_list_of_fractions_exp [x] = error ("RATIONALS_STEP_ADD_LIST_OF_FRACTIONS_EXP_EXCEPTION: Nothing to add")

  2630   | step_add_list_of_fractions_exp (xs)=

  2631     let

  2632         val den_list=termlist2denominators_exp (xs); (* list of denominators *)

  2633 	val (denom,var)=calc_lcm(den_list);      (* common denominator *)

  2634 	val den=factorize_den_exp(#1(den_list),denom,var); (* faktorisierter Nenner !!! *)

  2635     in

  2636 	com_den_exp(xs,denom,den,var)

  2637     end;

  2639 (* wird aktuell nicht mehr gebraucht, bei rückänderung schon

  2640 -------------------------------------------------------------

  2641 (* WN0210???SK brauch ma des überhaupt *)

  2642 fun step_add_list_of_fractions2 []=(Free("0",HOLogic.realT),[]:term list)

  2643   | step_add_list_of_fractions2 [x]=(x,[])

  2644   | step_add_list_of_fractions2 (xs) =

  2645     let

  2646         val den_list=termlist2denominators (xs); (* list of denominators *)

  2647 	val (denom,var)=calc_lcm(den_list);      (* common denominator *)

  2648 	val den=factorize_den(#1(den_list),denom,var);  (* faktorisierter Nenner !!! *)

  2649     in

  2650 	(

  2651 	 Const ("Rings.inverse_class.divide",[HOLogic.realT,HOLogic.realT]--->HOLogic.realT) $

  2652 	 com_den2(xs,denom, poly2term(denom,var)(*den*),var) $

  2653 	 poly2term(denom,var)

  2654 	,

  2655 	[]

  2656 	)

  2657     end;

  2659 (* WN0210???SK brauch ma des überhaupt *)

  2660 fun step_add_list_of_fractions2_exp []=(Free("0",HOLogic.realT),[]:term list)

  2661   | step_add_list_of_fractions2_exp [x]=(x,[])

  2662   | step_add_list_of_fractions2_exp (xs) =

  2663     let

  2664         val den_list=termlist2denominators_exp (xs); (* list of denominators *)

  2665 	val (denom,var)=calc_lcm(den_list);      (* common denominator *)

  2666 	val den=factorize_den_exp(#1(den_list),denom,var);  (* faktorisierter Nenner !!! *)

  2667     in

  2668 	(

  2669 	 Const ("Rings.inverse_class.divide",[HOLogic.realT,HOLogic.realT]--->HOLogic.realT) $

  2670 	 com_den_exp2(xs,denom, poly2term(denom,var)(*den*),var) $

  2671 	 poly2expanded(denom,var)

  2672 	,

  2673 	[]

  2674 	)

  2675     end;

  2676 ---------------------------------------------- *)

  2679 (*. converts a term, which contains severel terms seperated by +, into a list of these terms .*)

  2680 fun term2list (t as (Const("Rings.inverse_class.divide",_) $ _ $ _)) = [t]

  2681   | term2list (t as (Const("Atools.pow",_) $ _ $ _)) =

  2682     [Const("Rings.inverse_class.divide",[HOLogic.realT,HOLogic.realT]--->HOLogic.realT) $

  2683 	  t $ Free("1",HOLogic.realT)

  2684      ]

  2685   | term2list (t as (Free(_,_))) =

  2686     [Const("Rings.inverse_class.divide",[HOLogic.realT,HOLogic.realT]--->HOLogic.realT) $

  2687 	  t $  Free("1",HOLogic.realT)

  2688      ]

  2689   | term2list (t as (Const("Groups.times_class.times",_) $ _ $ _)) =

  2690     [Const("Rings.inverse_class.divide",[HOLogic.realT,HOLogic.realT]--->HOLogic.realT) $

  2691 	  t $ Free("1",HOLogic.realT)

  2692      ]

  2693   | term2list (Const("Groups.plus_class.plus",_) $ t1 $ t2) = term2list(t1) @ term2list(t2)

  2694   | term2list (Const("Groups.minus_class.minus",_) $ t1 $ t2) =

  2695     error ("RATIONALS_TERM2LIST_EXCEPTION: - not implemented yet")

  2696   | term2list _ = error ("RATIONALS_TERM2LIST_EXCEPTION: invalid term");

  2698 (*.factors out the gcd of nominator and denominator:

  2699    a/b = (a' * gcd)/(b' * gcd),  a,b,gcd  are poly[2].*)

  2700 fun factout_p_  (thy:theory) t = SOME (step_cancel t,[]:term list);

  2701 fun factout_ (thy:theory) t = SOME (step_cancel_expanded t,[]:term list);

  2703 (*.cancels a single fraction with normalform [2]

  2704    resulting in a canceled fraction [2], see factout_ .*)

  2705 fun cancel_p_ (thy:theory) t = (*WN.2.6.03 no rewrite -> NONE !*)

  2706     (let val (t',asm) = direct_cancel(*_expanded ... corrected MG.21.8.03*) t

  2707      in if t = t' then NONE else SOME (t',asm)

  2708      end) handle _ => NONE;

  2709 (*.the same as above with normalform [3]

  2710   val cancel_ :

  2711       theory ->        (*10.02 unused                                    *)

  2712       term -> 	       (*fraction in normalform [3]                      *)

  2713       (term * 	       (*fraction in normalform [3]                      *)

  2714        term list)      (*casual asumptions in normalform [3]             *)

  2715 	  option       (*NONE: the function is not applicable            *).*)

  2716 fun cancel_ (thy:theory) t = SOME (direct_cancel_expanded t) handle _ => NONE;

  2718 (*.transforms sums of at least 2 fractions [3] to

  2719    sums with the least common multiple as nominator.*)

  2720 fun common_nominator_p_ (thy:theory) t =

  2721 ((*tracing("### common_nominator_p_ called");*)

  2722     SOME (step_add_list_of_fractions(term2list(t))) handle _ => NONE

  2723 );

  2724 fun common_nominator_ (thy:theory) t =

  2725     SOME (step_add_list_of_fractions_exp(term2list(t))) handle _ => NONE;

  2727 (*.add 2 or more fractions

  2728 val add_fraction_p_ :

  2729       theory ->        (*10.02 unused                                    *)

  2730       term -> 	       (*2 or more fractions with normalform [2]         *)

  2731       (term * 	       (*one fraction with normalform [2]                *)

  2732        term list)      (*casual assumptions in normalform [2] WN0210???SK  *)

  2733 	  option       (*NONE: the function is not applicable            *).*)

  2734 fun add_fraction_p_ (thy:theory) t =

  2735 (tracing("### add_fraction_p_ called");

  2736     (let val ts = term2list t

  2737      in if 1 < length ts

  2738 	then SOME (add_list_of_fractions ts)

  2739 	else NONE (*error ("RATIONALS_ADD_EXCEPTION: nothing to add")*)

  2740      end) handle _ => NONE

  2741 );

  2742 (*.same as add_fraction_p_ but with normalform [3].*)

  2743 (*SOME (step_add_list_of_fractions2(term2list(t))); *)

  2744 fun add_fraction_ (thy:theory) t =

  2745     if length(term2list(t))>1

  2746     then SOME (add_list_of_fractions_exp(term2list(t))) handle _ => NONE

  2747     else (*error ("RATIONALS_ADD_FRACTION_EXCEPTION: nothing to add")*)

  2748 	NONE;

  2749 fun add_fraction_ (thy:theory) t =

  2750     (if 1 < length (term2list t)

  2751      then SOME (add_list_of_fractions_exp (term2list t))

  2752      else (*error ("RATIONALS_ADD_FRACTION_EXCEPTION: nothing to add")*)

  2753 	 NONE) handle _ => NONE;

  2755 (*SOME (step_add_list_of_fractions2_exp(term2list(t))); *)

  2757 (*. brings the term into a normal form .*)

  2758 fun norm_rational_ (thy:theory) t =

  2759     SOME (add_list_of_fractions(term2list(t))) handle _ => NONE;

  2760 fun norm_expanded_rat_ (thy:theory) t =

  2761     SOME (add_list_of_fractions_exp(term2list(t))) handle _ => NONE;

  2764 (*.evaluates conditions in calculate_Rational.*)

  2765 (*make local with FIXX@ME result:term *term list*)

  2766 val calc_rat_erls = prep_rls(

  2767   Rls {id = "calc_rat_erls", preconds = [], rew_ord = ("dummy_ord",dummy_ord),

  2768 	 erls = e_rls, srls = Erls, calc = [], (*asm_thm = [], *)

  2769 	 rules =

  2770 	 [Calc ("op =",eval_equal "#equal_"),

  2771 	  Calc ("Atools.is'_const",eval_const "#is_const_"),

  2772 	  Thm ("not_true",num_str @{thm not_true}),

  2773 	  Thm ("not_false",num_str @{thm not_false})

  2774 	  ],

  2775 	 scr = EmptyScr});

  2778 (*.simplifies expressions with numerals;

  2779    does NOT rearrange the term by AC-rewriting; thus terms with variables

  2780    need to have constants to be commuted together respectively.*)

  2781 val calculate_Rational = prep_rls(

  2782     merge_rls "calculate_Rational"

  2783 	(Rls {id = "divide", preconds = [], rew_ord = ("dummy_ord",dummy_ord),

  2784 	      erls = calc_rat_erls, srls = Erls, (*asm_thm = [],*)

  2785 	      calc = [],

  2786 	      rules =

  2787 	      [Calc ("Rings.inverse_class.divide"  ,eval_cancel "#divide_e"),

  2789 	       Thm ("minus_divide_left",

  2790 		    num_str (@{thm minus_divide_left} RS @{thm sym})),

  2791 	       (*SYM - ?x / ?y = - (?x / ?y)  may come from subst*)

  2793 	       Thm ("rat_add",num_str @{thm rat_add}),

  2794 	       (*"[| a is_const; b is_const; c is_const; d is_const |] ==> \

  2795 		 \"a / c + b / d = (a * d) / (c * d) + (b * c ) / (d * c)"*)

  2796 	       Thm ("rat_add1",num_str @{thm rat_add1}),

  2797 	       (*"[| a is_const; b is_const; c is_const |] ==> \

  2798 		 \"a / c + b / c = (a + b) / c"*)

  2799 	       Thm ("rat_add2",num_str @{thm rat_add2}),

  2800 	       (*"[| ?a is_const; ?b is_const; ?c is_const |] ==> \

  2801 		 \?a / ?c + ?b = (?a + ?b * ?c) / ?c"*)

  2802 	       Thm ("rat_add3",num_str @{thm rat_add3}),

  2803 	       (*"[| a is_const; b is_const; c is_const |] ==> \

  2804 		 \"a + b / c = (a * c) / c + b / c"\

  2805 		 \.... is_const to be omitted here FIXME*)

  2807 	       Thm ("rat_mult",num_str @{thm rat_mult}),

  2808 	       (*a / b * (c / d) = a * c / (b * d)*)

  2809 	       Thm ("times_divide_eq_right",num_str @{thm times_divide_eq_right}),

  2810 	       (*?x * (?y / ?z) = ?x * ?y / ?z*)

  2811 	       Thm ("times_divide_eq_left",num_str @{thm times_divide_eq_left}),

  2812 	       (*?y / ?z * ?x = ?y * ?x / ?z*)

  2814 	       Thm ("real_divide_divide1",num_str @{thm real_divide_divide1}),

  2815 	       (*"?y ~= 0 ==> ?u / ?v / (?y / ?z) = ?u / ?v * (?z / ?y)"*)

  2816 	       Thm ("divide_divide_eq_left",num_str @{thm divide_divide_eq_left}),

  2817 	       (*"?x / ?y / ?z = ?x / (?y * ?z)"*)

  2819 	       Thm ("rat_power", num_str @{thm rat_power}),

  2820 	       (*"(?a / ?b) ^^^ ?n = ?a ^^^ ?n / ?b ^^^ ?n"*)

  2822 	       Thm ("mult_cross",num_str @{thm mult_cross}),

  2823 	       (*"[| b ~= 0; d ~= 0 |] ==> (a / b = c / d) = (a * d = b * c)*)

  2824 	       Thm ("mult_cross1",num_str @{thm mult_cross1}),

  2825 	       (*"   b ~= 0            ==> (a / b = c    ) = (a     = b * c)*)

  2826 	       Thm ("mult_cross2",num_str @{thm mult_cross2})

  2827 	       (*"           d ~= 0    ==> (a     = c / d) = (a * d =     c)*)

  2828 	       ], scr = EmptyScr})

  2829 	calculate_Poly);

  2831 (*("is_expanded", ("Rational.is'_expanded", eval_is_expanded ""))*)

  2832 fun eval_is_expanded (thmid:string) _

  2833 		       (t as (Const("Rational.is'_expanded", _) $ arg)) thy =

  2834     if is_expanded arg

  2835     then SOME (mk_thmid thmid ""

  2836 			((Syntax.string_of_term (thy2ctxt thy)) arg) "",

  2837 	       Trueprop $ (mk_equality (t, HOLogic.true_const)))

  2838     else SOME (mk_thmid thmid ""

  2839 			((Syntax.string_of_term (thy2ctxt thy)) arg) "",

  2840 	       Trueprop $ (mk_equality (t, HOLogic.false_const)))

  2841   | eval_is_expanded _ _ _ _ = NONE;

  2843 val rational_erls =

  2844     merge_rls "rational_erls" calculate_Rational

  2845 	      (append_rls "is_expanded" Atools_erls

  2846 			  [Calc ("Rational.is'_expanded", eval_is_expanded "")

  2847 			   ]);

  2850 (*.3 'reverse-rewrite-sets' for symbolic computation on rationals:

  2851  =================================================================

  2852  A[2] 'cancel_p': .

  2853  A[3] 'cancel': .

  2854  B[2] 'common_nominator_p': transforms summands in a term [2]

  2855          to fractions with the (least) common multiple as nominator.

  2856  B[3] 'norm_rational': normalizes arbitrary algebraic terms (without

  2857          radicals and transzendental functions) to one canceled fraction,

  2858 	 nominator and denominator in polynomial form.

  2860 In order to meet isac's requirements for interactive and stepwise calculation,

  2861 each 'reverse-rewerite-set' consists of an initialization for the interpreter

  2862 state and of 4 functions, each of which employs rewriting as much as possible.

  2863 The signature of these functions are the same in each 'reverse-rewrite-set'

  2864 respectively.*)

  2866 (* ************************************************************************* *)

  2868 local(*. cancel_p

  2869 ------------------------

  2870 cancels a single fraction consisting of two (uni- or multivariate)

  2871 polynomials WN0609???SK[2] into another such a fraction; examples:

  2873 	   a^2 + -1*b^2         a + b

  2874         -------------------- = ---------

  2875 	a^2 + -2*a*b + b^2     a + -1*b

  2877         a^2    a

  2878         --- = ---

  2879          a     1

  2881 Remark: the reverse ruleset does _NOT_ work properly with other input !.*)

  2882 (*WN020824 wir werden "uberlegen, wie wir ungeeignete inputs zur"uckweisen*)

  2884 val {rules, rew_ord=(_,ro),...} =

  2885     rep_rls (assoc_rls "make_polynomial");

  2886 (*WN060829 ... make_deriv does not terminate with 1st expl above,

  2887            see rational.sml --- investigate rulesets for cancel_p ---*)

  2888 val {rules, rew_ord=(_,ro),...} =

  2889     rep_rls (assoc_rls "rev_rew_p");

  2891 (*.init_state = fn : term -> istate

  2892 initialzies the state of the script interpreter. The state is:

  2894 type rrlsstate =      (*state for reverse rewriting*)

  2895      (term *          (*the current formula*)

  2896       term *          (*the final term*)

  2897       rule list       (*'reverse rule list' (#)*)

  2898 	    list *    (*may be serveral, eg. in norm_rational*)

  2899       (rule *         (*Thm (+ Thm generated from Calc) resulting in ...*)

  2900        (term *        (*... rewrite with ...*)

  2901 	term list))   (*... assumptions*)

  2902 	  list);      (*derivation from given term to normalform

  2903 		       in reverse order with sym_thm;

  2904                        (#) could be extracted from here by (map #1)*).*)

  2905 (* val {rules, rew_ord=(_,ro),...} =

  2906        rep_rls (assoc_rls "rev_rew_p")        (*USE ALWAYS, SEE val cancel_p*);

  2907    val (thy, eval_rls, ro) =(thy, Atools_erls, ro) (*..val cancel_p*);

  2908    val t = t;

  2909    *)

  2910 fun init_state thy eval_rls ro t =

  2911     let val SOME (t',_) = factout_p_ thy t

  2912         val SOME (t'',asm) = cancel_p_ thy t

  2913         val der = reverse_deriv thy eval_rls rules ro NONE t'

  2914         val der = der @ [(Thm ("real_mult_div_cancel2",

  2915 			       num_str @{thm real_mult_div_cancel2}),

  2916 			  (t'',asm))]

  2917         val rs = (distinct_Thm o (map #1)) der

  2918 	val rs = filter_out (eq_Thms ["sym_real_add_zero_left",

  2919 				      "sym_real_mult_0",

  2920 				      "sym_real_mult_1"

  2921 				      (*..insufficient,eg.make_Polynomial*)])rs

  2922     in (t,t'',[rs(*here only _ONE_ to ease locate_rule*)],der) end;

  2924 (*.locate_rule = fn : rule list -> term -> rule

  2925 		      -> (rule * (term * term list) option) list.

  2926   checks a rule R for being a cancel-rule, and if it is,

  2927   then return the list of rules (+ the terms they are rewriting to)

  2928   which need to be applied before R should be applied.

  2929   precondition: the rule is applicable to the argument-term.

  2930 arguments:

  2931   rule list: the reverse rule list

  2932   -> term  : ... to which the rule shall be applied

  2933   -> rule  : ... to be applied to term

  2934 value:

  2935   -> (rule           : a rule rewriting to ...

  2936       * (term        : ... the resulting term ...

  2937          * term list): ... with the assumptions ( //#0).

  2938       ) list         : there may be several such rules;

  2939 		       the list is empty, if the rule has nothing to do

  2940 		       with cancelation.*)

  2941 (* val () = ();

  2942    *)

  2943 fun locate_rule thy eval_rls ro [rs] t r =

  2944     if (id_of_thm r) mem (map (id_of_thm)) rs

  2945     then let val ropt =

  2946 		 rewrite_ thy ro eval_rls true (thm_of_thm r) t;

  2947 	 in case ropt of

  2948 		SOME ta => [(r, ta)]

  2949 	      | NONE => (tracing("### locate_rule:  rewrite "^

  2950 				 (id_of_thm r)^" "^(term2str t)^" = NONE");

  2951 			 []) end

  2952     else (tracing("### locate_rule:  "^(id_of_thm r)^" not mem rrls");[])

  2953   | locate_rule _ _ _ _ _ _ =

  2954     error ("locate_rule: doesnt match rev-sets in istate");

  2956 (*.next_rule = fn : rule list -> term -> rule option

  2957   for a given term return the next rules to be done for cancelling.

  2958 arguments:

  2959   rule list     : the reverse rule list

  2960   term          : the term for which ...

  2961 value:

  2962   -> rule option: ... this rule is appropriate for cancellation;

  2963 		  there may be no such rule (if the term is canceled already.*)

  2964 (* val thy = thy;

  2965    val Rrls {rew_ord=(_,ro),...} = cancel;

  2966    val ([rs],t) = (rss,f);

  2967    next_rule thy eval_rls ro [rs] t;(*eval fun next_rule ... before!*)

  2969    val (thy, [rs]) = (thy, revsets);

  2970    val Rrls {rew_ord=(_,ro),...} = cancel;

  2971    nex [rs] t;

  2972    *)

  2973 fun next_rule thy eval_rls ro [rs] t =

  2974     let val der = make_deriv thy eval_rls rs ro NONE t;

  2975     in case der of

  2976 (* val (_,r,_)::_ = der;

  2977    *)

  2978 	   (_,r,_)::_ => SOME r

  2979 	 | _ => NONE

  2980     end

  2981   | next_rule _ _ _ _ _ =

  2982     error ("next_rule: doesnt match rev-sets in istate");

  2984 (*.val attach_form = f : rule list -> term -> term

  2985 			 -> (rule * (term * term list)) list

  2986   checks an input term TI, if it may belong to a current cancellation, by

  2987   trying to derive it from the given term TG.

  2988 arguments:

  2989   term   : TG, the last one in the cancellation agreed upon by user + math-eng

  2990   -> term: TI, the next one input by the user

  2991 value:

  2992   -> (rule           : the rule to be applied in order to reach TI

  2993       * (term        : ... obtained by applying the rule ...

  2994          * term list): ... and the respective assumptions.

  2995       ) list         : there may be several such rules;

  2996                        the list is empty, if the users term does not belong

  2997 		       to a cancellation of the term last agreed upon.*)

  2998 fun attach_form (_:rule list list) (_:term) (_:term) = (*still missing*)

  2999     []:(rule * (term * term list)) list;

  3001 in

  3003 val cancel_p =

  3004     Rrls {id = "cancel_p", prepat=[],

  3005 	  rew_ord=("ord_make_polynomial",

  3006 		   ord_make_polynomial false thy),

  3007 	  erls = rational_erls,

  3008 	  calc = [("PLUS"    ,("Groups.plus_class.plus"        ,eval_binop "#add_")),

  3009 		  ("TIMES"   ,("Groups.times_class.times"        ,eval_binop "#mult_")),

  3010 		  ("DIVIDE" ,("Rings.inverse_class.divide"  ,eval_cancel "#divide_e")),

  3011 		  ("POWER"  ,("Atools.pow"  ,eval_binop "#power_"))],

  3012 	  (*asm_thm=[("real_mult_div_cancel2","")],*)

  3013 	  scr=Rfuns {init_state  = init_state thy Atools_erls ro,

  3014 		     normal_form = cancel_p_ thy,

  3015 		     locate_rule = locate_rule thy Atools_erls ro,

  3016 		     next_rule   = next_rule thy Atools_erls ro,

  3017 		     attach_form = attach_form}}

  3018 end;(*local*)

  3020 local(*.ad (1) 'cancel'

  3021 ------------------------------

  3022 cancels a single fraction consisting of two (uni- or multivariate)

  3023 polynomials WN0609???SK[3] into another such a fraction; examples:

  3025 	   a^2 - b^2           a + b

  3026         -------------------- = ---------

  3027 	a^2 - 2*a*b + b^2      a - *b

  3029 Remark: the reverse ruleset does _NOT_ work properly with other input !.*)

  3030 (*WN 24.8.02: wir werden "uberlegen, wie wir ungeeignete inputs zur"uckweisen*)

  3032 (*

  3033 val SOME (Rls {rules=rules,rew_ord=(_,ro),...}) =

  3034     assoc'(!ruleset',"expand_binoms");

  3035 *)

  3036 val {rules=rules,rew_ord=(_,ro),...} =

  3037     rep_rls (assoc_rls "expand_binoms");

  3038 val thy = thy;

  3040 fun init_state thy eval_rls ro t =

  3041     let val SOME (t',_) = factout_ thy t;

  3042         val SOME (t'',asm) = cancel_ thy t;

  3043         val der = reverse_deriv thy eval_rls rules ro NONE t';

  3044         val der = der @ [(Thm ("real_mult_div_cancel2",

  3045 			       num_str @{thm real_mult_div_cancel2}),

  3046 			  (t'',asm))]

  3047         val rs = map #1 der;

  3048     in (t,t'',[rs],der) end;

  3050 fun locate_rule thy eval_rls ro [rs] t r =

  3051     if (id_of_thm r) mem (map (id_of_thm)) rs

  3052     then let val ropt =

  3053 		 rewrite_ thy ro eval_rls true (thm_of_thm r) t;

  3054 	 in case ropt of

  3055 		SOME ta => [(r, ta)]

  3056 	      | NONE => (tracing("### locate_rule:  rewrite "^

  3057 				 (id_of_thm r)^" "^(term2str t)^" = NONE");

  3058 			 []) end

  3059     else (tracing("### locate_rule:  "^(id_of_thm r)^" not mem rrls");[])

  3060   | locate_rule _ _ _ _ _ _ =

  3061     error ("locate_rule: doesnt match rev-sets in istate");

  3063 fun next_rule thy eval_rls ro [rs] t =

  3064     let val der = make_deriv thy eval_rls rs ro NONE t;

  3065     in case der of

  3066 (* val (_,r,_)::_ = der;

  3067    *)

  3068 	   (_,r,_)::_ => SOME r

  3069 	 | _ => NONE

  3070     end

  3071   | next_rule _ _ _ _ _ =

  3072     error ("next_rule: doesnt match rev-sets in istate");

  3074 fun attach_form (_:rule list list) (_:term) (_:term) = (*still missing*)

  3075     []:(rule * (term * term list)) list;

  3077 (* pat matched with the current term gives an environment

  3078    (or not: hen the Rrls not applied);

  3079    if pre1 and pre2 evaluate to HOLogic.true_const in this environment,

  3080    then the Rrls is applied. *)

  3081 val pat = parse_patt thy "?r/?s";

  3082 val pre1 = parse_patt thy "?r is_expanded";

  3083 val pre2 = parse_patt thy "?s is_expanded";

  3084 val prepat = [([pre1, pre2], pat)];

  3086 in

  3089 val cancel =

  3090     Rrls {id = "cancel", prepat=prepat,

  3091 	  rew_ord=("ord_make_polynomial",

  3092 		   ord_make_polynomial false thy),

  3093 	  erls = rational_erls,

  3094 	  calc = [("PLUS"    ,("Groups.plus_class.plus"        ,eval_binop "#add_")),

  3095 		  ("TIMES"   ,("Groups.times_class.times"        ,eval_binop "#mult_")),

  3096 		  ("DIVIDE" ,("Rings.inverse_class.divide"  ,eval_cancel "#divide_e")),

  3097 		  ("POWER"  ,("Atools.pow"  ,eval_binop "#power_"))],

  3098 	  scr=Rfuns {init_state  = init_state thy Atools_erls ro,

  3099 		     normal_form = cancel_ thy,

  3100 		     locate_rule = locate_rule thy Atools_erls ro,

  3101 		     next_rule   = next_rule thy Atools_erls ro,

  3102 		     attach_form = attach_form}}

  3103 end;(*local*)

  3105 local(*.ad [2] 'common_nominator_p'

  3106 ---------------------------------

  3107 FIXME Beschreibung .*)

  3110 val {rules=rules,rew_ord=(_,ro),...} =

  3111     rep_rls (assoc_rls "make_polynomial");

  3112 (*WN060829 ... make_deriv does not terminate with 1st expl above,

  3113            see rational.sml --- investigate rulesets for cancel_p ---*)

  3114 val {rules, rew_ord=(_,ro),...} =

  3115     rep_rls (assoc_rls "rev_rew_p");

  3116 val thy = thy;

  3119 (*.common_nominator_p_ = fn : theory -> term -> (term * term list) option

  3120   as defined above*)

  3122 (*.init_state = fn : term -> istate

  3123 initialzies the state of the interactive interpreter. The state is:

  3125 type rrlsstate =      (*state for reverse rewriting*)

  3126      (term *          (*the current formula*)

  3127       term *          (*the final term*)

  3128       rule list       (*'reverse rule list' (#)*)

  3129 	    list *    (*may be serveral, eg. in norm_rational*)

  3130       (rule *         (*Thm (+ Thm generated from Calc) resulting in ...*)

  3131        (term *        (*... rewrite with ...*)

  3132 	term list))   (*... assumptions*)

  3133 	  list);      (*derivation from given term to normalform

  3134 		       in reverse order with sym_thm;

  3135                        (#) could be extracted from here by (map #1)*).*)

  3136 fun init_state thy eval_rls ro t =

  3137     let val SOME (t',_) = common_nominator_p_ thy t;

  3138         val SOME (t'',asm) = add_fraction_p_ thy t;

  3139         val der = reverse_deriv thy eval_rls rules ro NONE t';

  3140         val der = der @ [(Thm ("real_mult_div_cancel2",

  3141 			       num_str @{thm real_mult_div_cancel2}),

  3142 			  (t'',asm))]

  3143         val rs = (distinct_Thm o (map #1)) der;

  3144 	val rs = filter_out (eq_Thms ["sym_real_add_zero_left",

  3145 				      "sym_real_mult_0",

  3146 				      "sym_real_mult_1"]) rs;

  3147     in (t,t'',[rs(*here only _ONE_*)],der) end;

  3149 (* use"knowledge/Rational.ML";

  3150    *)

  3152 (*.locate_rule = fn : rule list -> term -> rule

  3153 		      -> (rule * (term * term list) option) list.

  3154   checks a rule R for being a cancel-rule, and if it is,

  3155   then return the list of rules (+ the terms they are rewriting to)

  3156   which need to be applied before R should be applied.

  3157   precondition: the rule is applicable to the argument-term.

  3158 arguments:

  3159   rule list: the reverse rule list

  3160   -> term  : ... to which the rule shall be applied

  3161   -> rule  : ... to be applied to term

  3162 value:

  3163   -> (rule           : a rule rewriting to ...

  3164       * (term        : ... the resulting term ...

  3165          * term list): ... with the assumptions ( //#0).

  3166       ) list         : there may be several such rules;

  3167 		       the list is empty, if the rule has nothing to do

  3168 		       with cancelation.*)

  3169 (* val () = ();

  3170    *)

  3171 fun locate_rule thy eval_rls ro [rs] t r =

  3172     if (id_of_thm r) mem (map (id_of_thm)) rs

  3173     then let val ropt =

  3174 		 rewrite_ thy ro eval_rls true (thm_of_thm r) t;

  3175 	 in case ropt of

  3176 		SOME ta => [(r, ta)]

  3177 	      | NONE => (tracing("### locate_rule:  rewrite "^

  3178 				 (id_of_thm r)^" "^(term2str t)^" = NONE");

  3179 			 []) end

  3180     else (tracing("### locate_rule:  "^(id_of_thm r)^" not mem rrls");[])

  3181   | locate_rule _ _ _ _ _ _ =

  3182     error ("locate_rule: doesnt match rev-sets in istate");

  3184 (*.next_rule = fn : rule list -> term -> rule option

  3185   for a given term return the next rules to be done for cancelling.

  3186 arguments:

  3187   rule list     : the reverse rule list

  3188   term          : the term for which ...

  3189 value:

  3190   -> rule option: ... this rule is appropriate for cancellation;

  3191 		  there may be no such rule (if the term is canceled already.*)

  3192 (* val thy = thy;

  3193    val Rrls {rew_ord=(_,ro),...} = cancel;

  3194    val ([rs],t) = (rss,f);

  3195    next_rule thy eval_rls ro [rs] t;(*eval fun next_rule ... before!*)

  3197    val (thy, [rs]) = (thy, revsets);

  3198    val Rrls {rew_ord=(_,ro),...} = cancel;

  3199    nex [rs] t;

  3200    *)

  3201 fun next_rule thy eval_rls ro [rs] t =

  3202     let val der = make_deriv thy eval_rls rs ro NONE t;

  3203     in case der of

  3204 (* val (_,r,_)::_ = der;

  3205    *)

  3206 	   (_,r,_)::_ => SOME r

  3207 	 | _ => NONE

  3208     end

  3209   | next_rule _ _ _ _ _ =

  3210     error ("next_rule: doesnt match rev-sets in istate");

  3212 (*.val attach_form = f : rule list -> term -> term

  3213 			 -> (rule * (term * term list)) list

  3214   checks an input term TI, if it may belong to a current cancellation, by

  3215   trying to derive it from the given term TG.

  3216 arguments:

  3217   term   : TG, the last one in the cancellation agreed upon by user + math-eng

  3218   -> term: TI, the next one input by the user

  3219 value:

  3220   -> (rule           : the rule to be applied in order to reach TI

  3221       * (term        : ... obtained by applying the rule ...

  3222          * term list): ... and the respective assumptions.

  3223       ) list         : there may be several such rules;

  3224                        the list is empty, if the users term does not belong

  3225 		       to a cancellation of the term last agreed upon.*)

  3226 fun attach_form (_:rule list list) (_:term) (_:term) = (*still missing*)

  3227     []:(rule * (term * term list)) list;

  3229 (* if each pat* matches with the current term, the Rrls is applied

  3230    (there are no preconditions to be checked, they are HOLogic.true_const) *)

  3231 val pat0 = parse_patt thy "?r/?s+?u/?v";

  3232 val pat1 = parse_patt thy "?r/?s+?u   ";

  3233 val pat2 = parse_patt thy "?r   +?u/?v";

  3234 val prepat = [([HOLogic.true_const], pat0),

  3235 	      ([HOLogic.true_const], pat1),

  3236 	      ([HOLogic.true_const], pat2)];

  3237 in

  3239 (*11.02 schnelle L"osung f"ur RL: Bruch auch gek"urzt;

  3240   besser w"are: auf 1 gemeinsamen Bruchstrich, Nenner und Z"ahler unvereinfacht

  3241   dh. wie common_nominator_p_, aber auf 1 Bruchstrich*)

  3242 val common_nominator_p =

  3243     Rrls {id = "common_nominator_p", prepat=prepat,

  3244 	  rew_ord=("ord_make_polynomial",

  3245 		   ord_make_polynomial false thy),

  3246 	  erls = rational_erls,

  3247 	  calc = [("PLUS"    ,("Groups.plus_class.plus"        ,eval_binop "#add_")),

  3248 		  ("TIMES"   ,("Groups.times_class.times"        ,eval_binop "#mult_")),

  3249 		  ("DIVIDE" ,("Rings.inverse_class.divide"  ,eval_cancel "#divide_e")),

  3250 		  ("POWER"  ,("Atools.pow"  ,eval_binop "#power_"))],

  3251 	  scr=Rfuns {init_state  = init_state thy Atools_erls ro,

  3252 		     normal_form = add_fraction_p_ thy,(*FIXME.WN0211*)

  3253 		     locate_rule = locate_rule thy Atools_erls ro,

  3254 		     next_rule   = next_rule thy Atools_erls ro,

  3255 		     attach_form = attach_form}}

  3256 end;(*local*)

  3258 local(*.ad [2] 'common_nominator'

  3259 ---------------------------------

  3260 FIXME Beschreibung .*)

  3263 val {rules=rules,rew_ord=(_,ro),...} =

  3264     rep_rls (assoc_rls "make_polynomial");

  3265 val thy = thy;

  3268 (*.common_nominator_ = fn : theory -> term -> (term * term list) option

  3269   as defined above*)

  3271 (*.init_state = fn : term -> istate

  3272 initialzies the state of the interactive interpreter. The state is:

  3273 type rrlsstate =      (*state for reverse rewriting*)

  3274      (term *          (*the current formula*)

  3275       term *          (*the final term*)

  3276       rule list       (*'reverse rule list' (#)*)

  3277 	    list *    (*may be serveral, eg. in norm_rational*)

  3278       (rule *         (*Thm (+ Thm generated from Calc) resulting in ...*)

  3279        (term *        (*... rewrite with ...*)

  3280 	term list))   (*... assumptions*)

  3281 	  list);      (*derivation from given term to normalform

  3282 		       in reverse order with sym_thm;

  3283                        (#) could be extracted from here by (map #1)*).*)

  3284 fun init_state thy eval_rls ro t =

  3285     let val SOME (t',_) = common_nominator_ thy t;

  3286         val SOME (t'',asm) = add_fraction_ thy t;

  3287         val der = reverse_deriv thy eval_rls rules ro NONE t';

  3288         val der = der @ [(Thm ("real_mult_div_cancel2",

  3289 			       num_str @{thm real_mult_div_cancel2}),

  3290 			  (t'',asm))]

  3291         val rs = (distinct_Thm o (map #1)) der;

  3292 	val rs = filter_out (eq_Thms ["sym_real_add_zero_left",

  3293 				      "sym_real_mult_0",

  3294 				      "sym_real_mult_1"]) rs;

  3295     in (t,t'',[rs(*here only _ONE_*)],der) end;

  3297 (*.locate_rule = fn : rule list -> term -> rule

  3298 		      -> (rule * (term * term list) option) list.

  3299   checks a rule R for being a cancel-rule, and if it is,

  3300   then return the list of rules (+ the terms they are rewriting to)

  3301   which need to be applied before R should be applied.

  3302   precondition: the rule is applicable to the argument-term.

  3303 arguments:

  3304   rule list: the reverse rule list

  3305   -> term  : ... to which the rule shall be applied

  3306   -> rule  : ... to be applied to term

  3307 value:

  3308   -> (rule           : a rule rewriting to ...

  3309       * (term        : ... the resulting term ...

  3310          * term list): ... with the assumptions ( //#0).

  3311       ) list         : there may be several such rules;

  3312 		       the list is empty, if the rule has nothing to do

  3313 		       with cancelation.*)

  3314 (* val () = ();

  3315    *)

  3316 fun locate_rule thy eval_rls ro [rs] t r =

  3317     if (id_of_thm r) mem (map (id_of_thm)) rs

  3318     then let val ropt =

  3319 		 rewrite_ thy ro eval_rls true (thm_of_thm r) t;

  3320 	 in case ropt of

  3321 		SOME ta => [(r, ta)]

  3322 	      | NONE => (tracing("### locate_rule:  rewrite "^

  3323 				 (id_of_thm r)^" "^(term2str t)^" = NONE");

  3324 			 []) end

  3325     else (tracing("### locate_rule:  "^(id_of_thm r)^" not mem rrls");[])

  3326   | locate_rule _ _ _ _ _ _ =

  3327     error ("locate_rule: doesnt match rev-sets in istate");

  3329 (*.next_rule = fn : rule list -> term -> rule option

  3330   for a given term return the next rules to be done for cancelling.

  3331 arguments:

  3332   rule list     : the reverse rule list

  3333   term          : the term for which ...

  3334 value:

  3335   -> rule option: ... this rule is appropriate for cancellation;

  3336 		  there may be no such rule (if the term is canceled already.*)

  3337 (* val thy = thy;

  3338    val Rrls {rew_ord=(_,ro),...} = cancel;

  3339    val ([rs],t) = (rss,f);

  3340    next_rule thy eval_rls ro [rs] t;(*eval fun next_rule ... before!*)

  3342    val (thy, [rs]) = (thy, revsets);

  3343    val Rrls {rew_ord=(_,ro),...} = cancel_p;

  3344    nex [rs] t;

  3345    *)

  3346 fun next_rule thy eval_rls ro [rs] t =

  3347     let val der = make_deriv thy eval_rls rs ro NONE t;

  3348     in case der of

  3349 (* val (_,r,_)::_ = der;

  3350    *)

  3351 	   (_,r,_)::_ => SOME r

  3352 	 | _ => NONE

  3353     end

  3354   | next_rule _ _ _ _ _ =

  3355     error ("next_rule: doesnt match rev-sets in istate");

  3357 (*.val attach_form = f : rule list -> term -> term

  3358 			 -> (rule * (term * term list)) list

  3359   checks an input term TI, if it may belong to a current cancellation, by

  3360   trying to derive it from the given term TG.

  3361 arguments:

  3362   term   : TG, the last one in the cancellation agreed upon by user + math-eng

  3363   -> term: TI, the next one input by the user

  3364 value:

  3365   -> (rule           : the rule to be applied in order to reach TI

  3366       * (term        : ... obtained by applying the rule ...

  3367          * term list): ... and the respective assumptions.

  3368       ) list         : there may be several such rules;

  3369                        the list is empty, if the users term does not belong

  3370 		       to a cancellation of the term last agreed upon.*)

  3371 fun attach_form (_:rule list list) (_:term) (_:term) = (*still missing*)

  3372     []:(rule * (term * term list)) list;

  3374 (* if each pat* matches with the current term, the Rrls is applied

  3375    (there are no preconditions to be checked, they are HOLogic.true_const) *)

  3376 val pat0 =  parse_patt thy "?r/?s+?u/?v";

  3377 val pat01 = parse_patt thy "?r/?s-?u/?v";

  3378 val pat1 =  parse_patt thy "?r/?s+?u   ";

  3379 val pat11 = parse_patt thy "?r/?s-?u   ";

  3380 val pat2 =  parse_patt thy "?r   +?u/?v";

  3381 val pat21 = parse_patt thy "?r   -?u/?v";

  3382 val prepat = [([HOLogic.true_const], pat0),

  3383 	      ([HOLogic.true_const], pat01),

  3384 	      ([HOLogic.true_const], pat1),

  3385 	      ([HOLogic.true_const], pat11),

  3386 	      ([HOLogic.true_const], pat2),

  3387 	      ([HOLogic.true_const], pat21)];

  3390 in

  3392 val common_nominator =

  3393     Rrls {id = "common_nominator", prepat=prepat,

  3394 	  rew_ord=("ord_make_polynomial",

  3395 		   ord_make_polynomial false thy),

  3396 	  erls = rational_erls,

  3397 	  calc = [("PLUS"    ,("Groups.plus_class.plus"        ,eval_binop "#add_")),

  3398 		  ("TIMES"   ,("Groups.times_class.times"        ,eval_binop "#mult_")),

  3399 		  ("DIVIDE" ,("Rings.inverse_class.divide"  ,eval_cancel "#divide_e")),

  3400 		  ("POWER"  ,("Atools.pow"  ,eval_binop "#power_"))],

  3401 	  (*asm_thm=[("real_mult_div_cancel2","")],*)

  3402 	  scr=Rfuns {init_state  = init_state thy Atools_erls ro,

  3403 		     normal_form = add_fraction_ (*NOT common_nominator_*) thy,

  3404 		     locate_rule = locate_rule thy Atools_erls ro,

  3405 		     next_rule   = next_rule thy Atools_erls ro,

  3406 		     attach_form = attach_form}}

  3408 end;(*local*)

  3409 end;(*struct*)

  3411 open RationalI;

  3412 (*##*)

  3414 (*.the expression contains + - * ^ / only ?.*)

  3415 fun is_ratpolyexp (Free _) = true

  3416   | is_ratpolyexp (Const ("Groups.plus_class.plus",_) $ Free _ $ Free _) = true

  3417   | is_ratpolyexp (Const ("Groups.minus_class.minus",_) $ Free _ $ Free _) = true

  3418   | is_ratpolyexp (Const ("Groups.times_class.times",_) $ Free _ $ Free _) = true

  3419   | is_ratpolyexp (Const ("Atools.pow",_) $ Free _ $ Free _) = true

  3420   | is_ratpolyexp (Const ("Rings.inverse_class.divide",_) $ Free _ $ Free _) = true

  3421   | is_ratpolyexp (Const ("Groups.plus_class.plus",_) $ t1 $ t2) =

  3422                ((is_ratpolyexp t1) andalso (is_ratpolyexp t2))

  3423   | is_ratpolyexp (Const ("Groups.minus_class.minus",_) $ t1 $ t2) =

  3424                ((is_ratpolyexp t1) andalso (is_ratpolyexp t2))

  3425   | is_ratpolyexp (Const ("Groups.times_class.times",_) $ t1 $ t2) =

  3426                ((is_ratpolyexp t1) andalso (is_ratpolyexp t2))

  3427   | is_ratpolyexp (Const ("Atools.pow",_) $ t1 $ t2) =

  3428                ((is_ratpolyexp t1) andalso (is_ratpolyexp t2))

  3429   | is_ratpolyexp (Const ("Rings.inverse_class.divide",_) $ t1 $ t2) =

  3430                ((is_ratpolyexp t1) andalso (is_ratpolyexp t2))

  3431   | is_ratpolyexp _ = false;

  3433 (*("is_ratpolyexp", ("Rational.is'_ratpolyexp", eval_is_ratpolyexp ""))*)

  3434 fun eval_is_ratpolyexp (thmid:string) _

  3435 		       (t as (Const("Rational.is'_ratpolyexp", _) $ arg)) thy =

  3436     if is_ratpolyexp arg

  3437     then SOME (mk_thmid thmid ""

  3438 			((Syntax.string_of_term (thy2ctxt thy)) arg) "",

  3439 	       Trueprop $ (mk_equality (t, HOLogic.true_const)))

  3440     else SOME (mk_thmid thmid ""

  3441 			((Syntax.string_of_term (thy2ctxt thy)) arg) "",

  3442 	       Trueprop $ (mk_equality (t, HOLogic.false_const)))

  3443   | eval_is_ratpolyexp _ _ _ _ = NONE;

  3447 (*-------------------18.3.03 --> struct <-----------vvv--*)

  3448 val add_fractions_p = common_nominator_p; (*FIXXXME:eilig f"ur norm_Rational*)

  3450 (*WN100906 removed in favour of discard_minus = discard_minus_ formerly

  3451    discard binary minus, shift unary minus into -1*;

  3452    unary minus before numerals are put into the numeral by parsing;

  3453    contains absolute minimum of thms for context in norm_Rational

  3454 val discard_minus = prep_rls(

  3455   Rls {id = "discard_minus", preconds = [], rew_ord = ("dummy_ord",dummy_ord),

  3456       erls = e_rls, srls = Erls, calc = [],

  3457       rules = [Thm ("real_diff_minus", num_str @{thm real_diff_minus}),

  3458 	       (*"a - b = a + -1 * b"*)

  3459 	       Thm ("sym_real_mult_minus1",

  3460                      num_str (@{thm real_mult_minus1} RS @{thm sym}))

  3461 	       (*- ?z = "-1 * ?z"*)

  3462 	       ],

  3463       scr = EmptyScr

  3464       }):rls;

  3465  *)

  3467 (*erls for calculate_Rational; make local with FIXX@ME result:term *term list*)

  3468 val powers_erls = prep_rls(

  3469   Rls {id = "powers_erls", preconds = [], rew_ord = ("dummy_ord",dummy_ord),

  3470       erls = e_rls, srls = Erls, calc = [], (*asm_thm = [],*)

  3471       rules = [Calc ("Atools.is'_atom",eval_is_atom "#is_atom_"),

  3472 	       Calc ("Atools.is'_even",eval_is_even "#is_even_"),

  3473 	       Calc ("op <",eval_equ "#less_"),

  3474 	       Thm ("not_false", num_str @{thm not_false}),

  3475 	       Thm ("not_true", num_str @{thm not_true}),

  3476 	       Calc ("Groups.plus_class.plus",eval_binop "#add_")

  3477 	       ],

  3478       scr = EmptyScr

  3479       }:rls);

  3480 (*.all powers over + distributed; atoms over * collected, other distributed

  3481    contains absolute minimum of thms for context in norm_Rational .*)

  3482 val powers = prep_rls(

  3483   Rls {id = "powers", preconds = [], rew_ord = ("dummy_ord",dummy_ord),

  3484       erls = powers_erls, srls = Erls, calc = [], (*asm_thm = [],*)

  3485       rules = [Thm ("realpow_multI", num_str @{thm realpow_multI}),

  3486 	       (*"(r * s) ^^^ n = r ^^^ n * s ^^^ n"*)

  3487 	       Thm ("realpow_pow",num_str @{thm realpow_pow}),

  3488 	       (*"(a ^^^ b) ^^^ c = a ^^^ (b * c)"*)

  3489 	       Thm ("realpow_oneI",num_str @{thm realpow_oneI}),

  3490 	       (*"r ^^^ 1 = r"*)

  3491 	       Thm ("realpow_minus_even",num_str @{thm realpow_minus_even}),

  3492 	       (*"n is_even ==> (- r) ^^^ n = r ^^^ n" ?-->discard_minus?*)

  3493 	       Thm ("realpow_minus_odd",num_str @{thm realpow_minus_odd}),

  3494 	       (*"Not (n is_even) ==> (- r) ^^^ n = -1 * r ^^^ n"*)

  3496 	       (*----- collect atoms over * -----*)

  3497 	       Thm ("realpow_two_atom",num_str @{thm realpow_two_atom}),

  3498 	       (*"r is_atom ==> r * r = r ^^^ 2"*)

  3499 	       Thm ("realpow_plus_1",num_str @{thm realpow_plus_1}),

  3500 	       (*"r is_atom ==> r * r ^^^ n = r ^^^ (n + 1)"*)

  3501 	       Thm ("realpow_addI_atom",num_str @{thm realpow_addI_atom}),

  3502 	       (*"r is_atom ==> r ^^^ n * r ^^^ m = r ^^^ (n + m)"*)

  3504 	       (*----- distribute none-atoms -----*)

  3505 	       Thm ("realpow_def_atom",num_str @{thm realpow_def_atom}),

  3506 	       (*"[| 1 < n; not(r is_atom) |]==>r ^^^ n = r * r ^^^ (n + -1)"*)

  3507 	       Thm ("realpow_eq_oneI",num_str @{thm realpow_eq_oneI}),

  3508 	       (*"1 ^^^ n = 1"*)

  3509 	       Calc ("Groups.plus_class.plus",eval_binop "#add_")

  3510 	       ],

  3511       scr = EmptyScr

  3512       }:rls);

  3513 (*.contains absolute minimum of thms for context in norm_Rational.*)

  3514 val rat_mult_divide = prep_rls(

  3515   Rls {id = "rat_mult_divide", preconds = [],

  3516        rew_ord = ("dummy_ord",dummy_ord),

  3517       erls = e_rls, srls = Erls, calc = [], (*asm_thm = [],*)

  3518       rules = [Thm ("rat_mult",num_str @{thm rat_mult}),

  3519 	       (*(1)"?a / ?b * (?c / ?d) = ?a * ?c / (?b * ?d)"*)

  3520 	       Thm ("times_divide_eq_right",num_str @{thm times_divide_eq_right}),

  3521 	       (*(2)"?a * (?c / ?d) = ?a * ?c / ?d" must be [2],

  3522 	       otherwise inv.to a / b / c = ...*)

  3523 	       Thm ("times_divide_eq_left",num_str @{thm times_divide_eq_left}),

  3524 	       (*"?a / ?b * ?c = ?a * ?c / ?b" order weights x^^^n too much

  3525 		     and does not commute a / b * c ^^^ 2 !*)

  3527 	       Thm ("divide_divide_eq_right",

  3528                      num_str @{thm divide_divide_eq_right}),

  3529 	       (*"?x / (?y / ?z) = ?x * ?z / ?y"*)

  3530 	       Thm ("divide_divide_eq_left",

  3531                      num_str @{thm divide_divide_eq_left}),

  3532 	       (*"?x / ?y / ?z = ?x / (?y * ?z)"*)

  3533 	       Calc ("Rings.inverse_class.divide"  ,eval_cancel "#divide_e")

  3534 	       ],

  3535       scr = EmptyScr

  3536       }:rls);

  3538 (*.contains absolute minimum of thms for context in norm_Rational.*)

  3539 val reduce_0_1_2 = prep_rls(

  3540   Rls{id = "reduce_0_1_2", preconds = [], rew_ord = ("dummy_ord", dummy_ord),

  3541       erls = e_rls,srls = Erls,calc = [],(*asm_thm = [],*)

  3542       rules = [(*Thm ("divide_1",num_str @{thm divide_1}),

  3543 		 "?x / 1 = ?x" unnecess.for normalform*)

  3544 	       Thm ("mult_1_left",num_str @{thm mult_1_left}),

  3545 	       (*"1 * z = z"*)

  3546 	       (*Thm ("real_mult_minus1",num_str @{thm real_mult_minus1}),

  3547 	       "-1 * z = - z"*)

  3548 	       (*Thm ("real_minus_mult_cancel",num_str @{thm real_minus_mult_cancel}),

  3549 	       "- ?x * - ?y = ?x * ?y"*)

  3551 	       Thm ("mult_zero_left",num_str @{thm mult_zero_left}),

  3552 	       (*"0 * z = 0"*)

  3553 	       Thm ("add_0_left",num_str @{thm add_0_left}),

  3554 	       (*"0 + z = z"*)

  3555 	       (*Thm ("right_minus",num_str @{thm right_minus}),

  3556 	       "?z + - ?z = 0"*)

  3558 	       Thm ("sym_real_mult_2",

  3559                      num_str (@{thm real_mult_2} RS @{thm sym})),

  3560 	       (*"z1 + z1 = 2 * z1"*)

  3561 	       Thm ("real_mult_2_assoc",num_str @{thm real_mult_2_assoc}),

  3562 	       (*"z1 + (z1 + k) = 2 * z1 + k"*)

  3564 	       Thm ("divide_zero_left",num_str @{thm divide_zero_left})

  3565 	       (*"0 / ?x = 0"*)

  3566 	       ], scr = EmptyScr}:rls);

  3568 (*erls for calculate_Rational;

  3569   make local with FIXX@ME result:term *term list WN0609???SKMG*)

  3570 val norm_rat_erls = prep_rls(

  3571   Rls {id = "norm_rat_erls", preconds = [], rew_ord = ("dummy_ord",dummy_ord),

  3572       erls = e_rls, srls = Erls, calc = [], (*asm_thm = [],*)

  3573       rules = [Calc ("Atools.is'_const",eval_const "#is_const_")

  3574 	       ], scr = EmptyScr}:rls);

  3576 (*.consists of rls containing the absolute minimum of thms.*)

  3577 (*040209: this version has been used by RL for his equations,

  3578 which is now replaced by MGs version below

  3579 vvv OLD VERSION !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!*)

  3580 val norm_Rational = prep_rls(

  3581   Rls {id = "norm_Rational", preconds = [], rew_ord = ("dummy_ord",dummy_ord),

  3582       erls = norm_rat_erls, srls = Erls, calc = [],

  3583       rules = [(*sequence given by operator precedence*)

  3584 	       Rls_ discard_minus,

  3585 	       Rls_ powers,

  3586 	       Rls_ rat_mult_divide,

  3587 	       Rls_ expand,

  3588 	       Rls_ reduce_0_1_2,

  3589 	       (*^^^^^^^^^ from RL -- not the latest one vvvvvvvvv*)

  3590 	       Rls_ order_add_mult,

  3591 	       Rls_ collect_numerals,

  3592 	       Rls_ add_fractions_p,

  3593 	       Rls_ cancel_p

  3594 	       ],

  3595       scr = EmptyScr}:rls);

  3597 val norm_Rational_parenthesized = prep_rls(

  3598   Seq {id = "norm_Rational_parenthesized", preconds = []:term list,

  3599        rew_ord = ("dummy_ord", dummy_ord),

  3600       erls = Atools_erls, srls = Erls,

  3601       calc = [],

  3602       rules = [Rls_  norm_Rational, (*from RL -- not the latest one*)

  3603 	       Rls_ discard_parentheses

  3604 	       ],

  3605       scr = EmptyScr}:rls);

  3607 (*WN030318???SK: simplifies all but cancel and common_nominator*)

  3608 val simplify_rational =

  3609     merge_rls "simplify_rational" expand_binoms

  3610     (append_rls "divide" calculate_Rational

  3611 		[Thm ("divide_1",num_str @{thm divide_1}),

  3612 		 (*"?x / 1 = ?x"*)

  3613 		 Thm ("rat_mult",num_str @{thm rat_mult}),

  3614 		 (*(1)"?a / ?b * (?c / ?d) = ?a * ?c / (?b * ?d)"*)

  3615 		 Thm ("times_divide_eq_right",num_str @{thm times_divide_eq_right}),

  3616 		 (*(2)"?a * (?c / ?d) = ?a * ?c / ?d" must be [2],

  3617 		 otherwise inv.to a / b / c = ...*)

  3618 		 Thm ("times_divide_eq_left",num_str @{thm times_divide_eq_left}),

  3619 		 (*"?a / ?b * ?c = ?a * ?c / ?b"*)

  3620 		 Thm ("add_minus",num_str @{thm add_minus}),

  3621 		 (*"?a + ?b - ?b = ?a"*)

  3622 		 Thm ("add_minus1",num_str @{thm add_minus1}),

  3623 		 (*"?a - ?b + ?b = ?a"*)

  3624 		 Thm ("divide_minus1",num_str @{thm divide_minus1})

  3625 		 (*"?x / -1 = - ?x"*)

  3626 (*

  3627 ,

  3628 		 Thm ("",num_str @{thm })

  3629 *)

  3630 		 ]);

  3632 (*---------vvv-------------MG ab 1.07.2003--------------vvv-----------*)

  3634 (* ------------------------------------------------------------------ *)

  3635 (*                  Simplifier für beliebige Buchterme                *)

  3636 (* ------------------------------------------------------------------ *)

  3637 (*----------------------- norm_Rational_mg ---------------------------*)

  3638 (*. description of the simplifier see MG-DA.p.56ff .*)

  3639 (* ------------------------------------------------------------------- *)

  3640 val common_nominator_p_rls = prep_rls(

  3641   Rls {id = "common_nominator_p_rls", preconds = [],

  3642        rew_ord = ("dummy_ord",dummy_ord),

  3643 	 erls = e_rls, srls = Erls, calc = [],

  3644 	 rules =

  3645 	 [Rls_ common_nominator_p

  3646 	  (*FIXME.WN0401 ? redesign Rrls - use exhaustively on a term ?

  3647 	    FIXME.WN0510 unnecessary nesting: introduce RRls_ : rls -> rule*)

  3648 	  ],

  3649 	 scr = EmptyScr});

  3650 (* ------------------------------------------------------------------- *)

  3651 val cancel_p_rls = prep_rls(

  3652   Rls {id = "cancel_p_rls", preconds = [],

  3653        rew_ord = ("dummy_ord",dummy_ord),

  3654 	 erls = e_rls, srls = Erls, calc = [],

  3655 	 rules =

  3656 	 [Rls_ cancel_p

  3657 	  (*FIXME.WN.0401 ? redesign Rrls - use exhaustively on a term ?*)

  3658 	  ],

  3659 	 scr = EmptyScr});

  3660 (* -------------------------------------------------------------------- *)

  3661 (*. makes 'normal' fractions; 'is_polyexp' inhibits double fractions;

  3662     used in initial part norm_Rational_mg, see example DA-M02-main.p.60.*)

  3663 val rat_mult_poly = prep_rls(

  3664   Rls {id = "rat_mult_poly", preconds = [],

  3665        rew_ord = ("dummy_ord",dummy_ord),

  3666 	 erls =  append_rls "e_rls-is_polyexp" e_rls

  3667 	         [Calc ("Poly.is'_polyexp", eval_is_polyexp "")],

  3668 	 srls = Erls, calc = [],

  3669 	 rules =

  3670 	 [Thm ("rat_mult_poly_l",num_str @{thm rat_mult_poly_l}),

  3671 	  (*"?c is_polyexp ==> ?c * (?a / ?b) = ?c * ?a / ?b"*)

  3672 	  Thm ("rat_mult_poly_r",num_str @{thm rat_mult_poly_r})

  3673 	  (*"?c is_polyexp ==> ?a / ?b * ?c = ?a * ?c / ?b"*)

  3674 	  ],

  3675 	 scr = EmptyScr});

  3677 (* ------------------------------------------------------------------ *)

  3678 (*. makes 'normal' fractions; 'is_polyexp' inhibits double fractions;

  3679     used in looping part norm_Rational_rls, see example DA-M02-main.p.60

  3680     .. WHERE THE LATTER DOES ALWAYS WORK, BECAUSE erls = e_rls,

  3681     I.E. THE RESPECTIVE ASSUMPTION IS STORED AND Thm APPLIED; WN051028

  3682     ... WN0609???MG.*)

  3683 val rat_mult_div_pow = prep_rls(

  3684   Rls {id = "rat_mult_div_pow", preconds = [],

  3685        rew_ord = ("dummy_ord",dummy_ord),

  3686        erls = e_rls,

  3687        (*FIXME.WN051028 append_rls "e_rls-is_polyexp" e_rls

  3688 			[Calc ("Poly.is'_polyexp", eval_is_polyexp "")],

  3689          with this correction ^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^ we get

  3690 	 error "rational.sml.sml: diff.behav. in norm_Rational_mg 29" etc.

  3691          thus we decided to go on with this flaw*)

  3692 		 srls = Erls, calc = [],

  3693       rules = [Thm ("rat_mult",num_str @{thm rat_mult}),

  3694 	       (*"?a / ?b * (?c / ?d) = ?a * ?c / (?b * ?d)"*)

  3695 	       Thm ("rat_mult_poly_l",num_str @{thm rat_mult_poly_l}),

  3696 	       (*"?c is_polyexp ==> ?c * (?a / ?b) = ?c * ?a / ?b"*)

  3697 	       Thm ("rat_mult_poly_r",num_str @{thm rat_mult_poly_r}),

  3698 	       (*"?c is_polyexp ==> ?a / ?b * ?c = ?a * ?c / ?b"*)

  3700 	       Thm ("real_divide_divide1_mg",

  3701                      num_str @{thm real_divide_divide1_mg}),

  3702 	       (*"y ~= 0 ==> (u / v) / (y / z) = (u * z) / (y * v)"*)

  3703 	       Thm ("divide_divide_eq_right",

  3704                      num_str @{thm divide_divide_eq_right}),

  3705 	       (*"?x / (?y / ?z) = ?x * ?z / ?y"*)

  3706 	       Thm ("divide_divide_eq_left",

  3707                      num_str @{thm divide_divide_eq_left}),

  3708 	       (*"?x / ?y / ?z = ?x / (?y * ?z)"*)

  3709 	       Calc ("Rings.inverse_class.divide"  ,eval_cancel "#divide_e"),

  3711 	       Thm ("rat_power", num_str @{thm rat_power})

  3712 		(*"(?a / ?b) ^^^ ?n = ?a ^^^ ?n / ?b ^^^ ?n"*)

  3713 	       ],

  3714       scr = EmptyScr}:rls);

  3715 (* ------------------------------------------------------------------ *)

  3716 val rat_reduce_1 = prep_rls(

  3717   Rls {id = "rat_reduce_1", preconds = [],

  3718        rew_ord = ("dummy_ord",dummy_ord),

  3719        erls = e_rls, srls = Erls, calc = [],

  3720        rules = [Thm ("divide_1",num_str @{thm divide_1}),

  3721 		(*"?x / 1 = ?x"*)

  3722 		Thm ("mult_1_left",num_str @{thm mult_1_left})

  3723 		(*"1 * z = z"*)

  3724 		],

  3725        scr = EmptyScr}:rls);

  3726 (* ------------------------------------------------------------------ *)

  3727 (*. looping part of norm_Rational(*_mg*) .*)

  3728 val norm_Rational_rls = prep_rls(

  3729    Rls {id = "norm_Rational_rls", preconds = [],

  3730        rew_ord = ("dummy_ord",dummy_ord),

  3731        erls = norm_rat_erls, srls = Erls, calc = [],

  3732        rules = [Rls_ common_nominator_p_rls,

  3733 		Rls_ rat_mult_div_pow,

  3734 		Rls_ make_rat_poly_with_parentheses,

  3735 		Rls_ cancel_p_rls,(*FIXME:cancel_p does NOT order sometimes*)

  3736 		Rls_ rat_reduce_1

  3737 		],

  3738        scr = EmptyScr}:rls);

  3739 (* ------------------------------------------------------------------ *)

  3740 (*040109 'norm_Rational'(by RL) replaced by 'norm_Rational_mg'(MG)

  3741  just be renaming:*)

  3742 val norm_Rational(*_mg*) = prep_rls(

  3743    Seq {id = "norm_Rational"(*_mg*), preconds = [],

  3744        rew_ord = ("dummy_ord",dummy_ord),

  3745        erls = norm_rat_erls, srls = Erls, calc = [],

  3746        rules = [Rls_ discard_minus,

  3747 		Rls_ rat_mult_poly,(* removes double fractions like a/b/c    *)

  3748 		Rls_ make_rat_poly_with_parentheses, (*WN0510 also in(#)below*)

  3749 		Rls_ cancel_p_rls, (*FIXME.MG:cancel_p does NOT order sometim*)

  3750 		Rls_ norm_Rational_rls,   (* the main rls, looping (#)       *)

  3751 		Rls_ discard_parentheses1 (* mult only                       *)

  3752 		],

  3753        scr = EmptyScr}:rls);

  3754 (* ------------------------------------------------------------------ *)

  3756 ruleset' := overwritelthy @{theory} (!ruleset',

  3757   [("calculate_Rational", calculate_Rational),

  3758    ("calc_rat_erls",calc_rat_erls),

  3759    ("rational_erls", rational_erls),

  3760    ("cancel_p", cancel_p),

  3761    ("cancel", cancel),

  3762    ("common_nominator_p", common_nominator_p),

  3763    ("common_nominator_p_rls", common_nominator_p_rls),

  3764    ("common_nominator"  , common_nominator),

  3765    ("discard_minus", discard_minus),

  3766    ("powers_erls", powers_erls),

  3767    ("powers", powers),

  3768    ("rat_mult_divide", rat_mult_divide),

  3769    ("reduce_0_1_2", reduce_0_1_2),

  3770    ("rat_reduce_1", rat_reduce_1),

  3771    ("norm_rat_erls", norm_rat_erls),

  3772    ("norm_Rational", norm_Rational),

  3773    ("norm_Rational_rls", norm_Rational_rls),

  3774    ("norm_Rational_parenthesized", norm_Rational_parenthesized),

  3775    ("rat_mult_poly", rat_mult_poly),

  3776    ("rat_mult_div_pow", rat_mult_div_pow),

  3777    ("cancel_p_rls", cancel_p_rls)

  3778    ]);

  3780 calclist':= overwritel (!calclist',

  3781    [("is_expanded", ("Rational.is'_expanded", eval_is_expanded ""))

  3782     ]);

  3784 (** problems **)

  3786 store_pbt

  3787  (prep_pbt thy "pbl_simp_rat" [] e_pblID

  3788  (["rational","simplification"],

  3789   [("#Given" ,["TERM t_t"]),

  3790    ("#Where" ,["t_t is_ratpolyexp"]),

  3791    ("#Find"  ,["normalform n_n"])

  3792   ],

  3793   append_rls "e_rls" e_rls [(*for preds in where_*)],

  3794   SOME "Simplifyt_t",

  3795   [["simplification","of_rationals"]]));

  3797 (** methods **)

  3798 val --------------------------------------------------+++++-+++++ = "00000";

  3800 (*WN061025 this methods script is copied from (auto-generated) script

  3801   of norm_Rational in order to ease repair on inform*)

  3802 store_met

  3803     (prep_met thy "met_simp_rat" [] e_metID

  3804 	      (["simplification","of_rationals"],

  3805 	       [("#Given" ,["TERM t_t"]),

  3806 		("#Where" ,["t_t is_ratpolyexp"]),

  3807 		("#Find"  ,["normalform n_n"])

  3808 		],

  3809 	       {rew_ord'="tless_true",

  3810 		rls' = e_rls,

  3811 		calc = [], srls = e_rls,

  3812 		prls = append_rls "simplification_of_rationals_prls" e_rls

  3813 				[(*for preds in where_*)

  3814 				 Calc ("Rational.is'_ratpolyexp",

  3815 				       eval_is_ratpolyexp "")],

  3816 		crls = e_rls, nrls = norm_Rational_rls},

  3817 "Script SimplifyScript (t_t::real) =                             " ^

  3818 "  ((Try (Rewrite_Set discard_minus False) @@                   " ^

  3819 "    Try (Rewrite_Set rat_mult_poly False) @@                    " ^

  3820 "    Try (Rewrite_Set make_rat_poly_with_parentheses False) @@   " ^

  3821 "    Try (Rewrite_Set cancel_p_rls False) @@                     " ^

  3822 "    (Repeat                                                     " ^

  3823 "     ((Try (Rewrite_Set common_nominator_p_rls False) @@        " ^

  3824 "       Try (Rewrite_Set rat_mult_div_pow False) @@              " ^

  3825 "       Try (Rewrite_Set make_rat_poly_with_parentheses False) @@" ^

  3826 "       Try (Rewrite_Set cancel_p_rls False) @@                  " ^

  3827 "       Try (Rewrite_Set rat_reduce_1 False)))) @@               " ^

  3828 "    Try (Rewrite_Set discard_parentheses1 False))               " ^

  3829 "   t_t)"

  3830 	       ));

  3832 *}

  3834 end

author	Walther Neuper <neuper@ist.tugraz.at>
	Fri, 01 Oct 2010 10:23:38 +0200
branch	isac-update-Isa09-2
changeset 38036	02a9b2540eb7
parent 38034	928cebc9c4aa
child 38037	a51a70334191
permissions	-rw-r--r--