1.1 --- a/src/Tools/isac/IsacKnowledge/Rational.ML Wed Aug 25 15:15:01 2010 +0200
1.2 +++ /dev/null Thu Jan 01 00:00:00 1970 +0000
1.3 @@ -1,3786 +0,0 @@
1.4 -(*.calculate in rationals: gcd, lcm, etc.
1.5 - (c) Stefan Karnel 2002
1.6 - Institute for Mathematics D and Institute for Software Technology,
1.7 - TU-Graz SS 2002
1.8 - Use is subject to license terms.
1.9 -
1.10 -use"IsacKnowledge/Rational.ML";
1.11 -use"Rational.ML";
1.12 -
1.13 -remove_thy"Rational";
1.14 -use_thy"IsacKnowledge/Isac";
1.15 -****************************************************************.*)
1.16 -
1.17 -(*.*****************************************************************
1.18 - Remark on notions in the documentation below:
1.19 - referring to the remark on 'polynomials' in Poly.sml we use
1.20 - [2] 'polynomial' normalform (Polynom)
1.21 - [3] 'expanded_term' normalform (Ausmultiplizierter Term),
1.22 - where normalform [2] is a special case of [3], i.e. [3] implies [2].
1.23 - Instead of
1.24 - 'fraction with numerator and nominator both in normalform [2]'
1.25 - 'fraction with numerator and nominator both in normalform [3]'
1.26 - we say:
1.27 - 'fraction in normalform [2]'
1.28 - 'fraction in normalform [3]'
1.29 - or
1.30 - 'fraction [2]'
1.31 - 'fraction [3]'.
1.32 - a 'simple fraction' is a term with '/' as outmost operator and
1.33 - numerator and nominator in normalform [2] or [3].
1.34 -****************************************************************.*)
1.35 -
1.36 -signature RATIONALI =
1.37 -sig
1.38 - type mv_monom
1.39 - type mv_poly
1.40 - val add_fraction_ : theory -> term -> (term * term list) option
1.41 - val add_fraction_p_ : theory -> term -> (term * term list) option
1.42 - val calculate_Rational : rls
1.43 - val calc_rat_erls:rls
1.44 - val cancel : rls
1.45 - val cancel_ : theory -> term -> (term * term list) option
1.46 - val cancel_p : rls
1.47 - val cancel_p_ : theory -> term -> (term * term list) option
1.48 - val common_nominator : rls
1.49 - val common_nominator_ : theory -> term -> (term * term list) option
1.50 - val common_nominator_p : rls
1.51 - val common_nominator_p_ : theory -> term -> (term * term list) option
1.52 - val eval_is_expanded : string -> 'a -> term -> theory ->
1.53 - (string * term) option
1.54 - val expanded2polynomial : term -> term option
1.55 - val factout_ : theory -> term -> (term * term list) option
1.56 - val factout_p_ : theory -> term -> (term * term list) option
1.57 - val is_expanded : term -> bool
1.58 - val is_polynomial : term -> bool
1.59 -
1.60 - val mv_gcd : (int * int list) list -> mv_poly -> mv_poly
1.61 - val mv_lcm : mv_poly -> mv_poly -> mv_poly
1.62 -
1.63 - val norm_expanded_rat_ : theory -> term -> (term * term list) option
1.64 -(*WN0602.2.6.pull into struct !!!
1.65 - val norm_Rational : rls(*.normalizes an arbitrary rational term without
1.66 - roots into a simple and canceled fraction
1.67 - with normalform [2].*)
1.68 -*)
1.69 -(*val norm_rational_p : 19.10.02 missing FIXXXXXXXXXXXXME
1.70 - rls (*.normalizes an rational term [2] without
1.71 - roots into a simple and canceled fraction
1.72 - with normalform [2].*)
1.73 -*)
1.74 - val norm_rational_ : theory -> term -> (term * term list) option
1.75 - val polynomial2expanded : term -> term option
1.76 - val rational_erls :
1.77 - rls (*.evaluates an arbitrary rational term with numerals.*)
1.78 -
1.79 -(*WN0210???SK: fehlen Funktionen, die exportiert werden sollen ? *)
1.80 -end
1.81 -
1.82 -(*.**************************************************************************
1.83 -survey on the functions
1.84 -~~~~~~~~~~~~~~~~~~~~~~~
1.85 - [2] 'polynomial' :rls | [3]'expanded_term':rls
1.86 ---------------------:------------------+-------------------:-----------------
1.87 - factout_p_ : | factout_ :
1.88 - cancel_p_ : | cancel_ :
1.89 - :cancel_p | :cancel
1.90 ---------------------:------------------+-------------------:-----------------
1.91 - common_nominator_p_: | common_nominator_ :
1.92 - :common_nominator_p| :common_nominator
1.93 - add_fraction_p_ : | add_fraction_ :
1.94 ---------------------:------------------+-------------------:-----------------
1.95 -???SK :norm_rational_p | :norm_rational
1.96 -
1.97 -This survey shows only the principal functions for reuse, and the identifiers
1.98 -of the rls exported. The list below shows some more useful functions.
1.99 -
1.100 -
1.101 -conversion from Isabelle-term to internal representation
1.102 -~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
1.103 -
1.104 -... BITTE FORTSETZEN ...
1.105 -
1.106 -polynomial2expanded = ...
1.107 -expanded2polynomial = ...
1.108 -
1.109 -remark: polynomial2expanded o expanded2polynomial = I,
1.110 - where 'o' is function chaining, and 'I' is identity WN0210???SK
1.111 -
1.112 -functions for greatest common divisor and canceling
1.113 -~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
1.114 -mv_gcd
1.115 -factout_
1.116 -factout_p_
1.117 -cancel_
1.118 -cancel_p_
1.119 -
1.120 -functions for least common multiple and addition of fractions
1.121 -~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
1.122 -mv_lcm
1.123 -common_nominator_
1.124 -common_nominator_p_
1.125 -add_fraction_ (*.add 2 or more fractions.*)
1.126 -add_fraction_p_ (*.add 2 or more fractions.*)
1.127 -
1.128 -functions for normalform of rationals
1.129 -~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
1.130 -WN0210???SK interne Funktionen f"ur norm_rational:
1.131 - schaffen diese SML-Funktionen wirklich ganz allgemeine Terme ?
1.132 -
1.133 -norm_rational_
1.134 -norm_expanded_rat_
1.135 -
1.136 -**************************************************************************.*)
1.137 -
1.138 -
1.139 -(*##*)
1.140 -structure RationalI : RATIONALI =
1.141 -struct
1.142 -(*##*)
1.143 -
1.144 -infix mem ins union; (*WN100819 updating to Isabelle2009-2*)
1.145 -fun x mem [] = false
1.146 - | x mem (y :: ys) = x = y orelse x mem ys;
1.147 -fun (x ins xs) = if x mem xs then xs else x :: xs;
1.148 -fun xs union [] = xs
1.149 - | [] union ys = ys
1.150 - | (x :: xs) union ys = xs union (x ins ys);
1.151 -
1.152 -(*. gcd of integers .*)
1.153 -(* die gcd Funktion von Isabelle funktioniert nicht richtig !!! *)
1.154 -fun gcd_int a b = if b=0 then a
1.155 - else gcd_int b (a mod b);
1.156 -
1.157 -(*. univariate polynomials (uv) .*)
1.158 -(*. univariate polynomials are represented as a list of the coefficent in reverse maximum degree order .*)
1.159 -(*. 5 * x^5 + 4 * x^3 + 2 * x^2 + x + 19 => [19,1,2,4,0,5] .*)
1.160 -type uv_poly = int list;
1.161 -
1.162 -(*. adds two uv polynomials .*)
1.163 -fun uv_mod_add_poly ([]:uv_poly,p2:uv_poly) = p2:uv_poly
1.164 - | uv_mod_add_poly (p1,[]) = p1
1.165 - | uv_mod_add_poly (x::p1,y::p2) = (x+y)::(uv_mod_add_poly(p1,p2));
1.166 -
1.167 -(*. multiplies a uv polynomial with a skalar s .*)
1.168 -fun uv_mod_smul_poly ([]:uv_poly,s:int) = []:uv_poly
1.169 - | uv_mod_smul_poly (x::p,s) = (x*s)::(uv_mod_smul_poly(p,s));
1.170 -
1.171 -(*. calculates the remainder of a polynomial divided by a skalar s .*)
1.172 -fun uv_mod_rem_poly ([]:uv_poly,s) = []:uv_poly
1.173 - | uv_mod_rem_poly (x::p,s) = (x mod s)::(uv_mod_smul_poly(p,s));
1.174 -
1.175 -(*. calculates the degree of a uv polynomial .*)
1.176 -fun uv_mod_deg ([]:uv_poly) = 0
1.177 - | uv_mod_deg p = length(p)-1;
1.178 -
1.179 -(*. calculates the remainder of x/p and represents it as value between -p/2 and p/2 .*)
1.180 -fun uv_mod_mod2(x,p)=
1.181 - let
1.182 - val y=(x mod p);
1.183 - in
1.184 - if (y)>(p div 2) then (y)-p else
1.185 - (
1.186 - if (y)<(~p div 2) then p+(y) else (y)
1.187 - )
1.188 - end;
1.189 -
1.190 -(*.calculates the remainder for each element of a integer list divided by p.*)
1.191 -fun uv_mod_list_modp [] p = []
1.192 - | uv_mod_list_modp (x::xs) p = (uv_mod_mod2(x,p))::(uv_mod_list_modp xs p);
1.193 -
1.194 -(*. appends an integer at the end of a integer list .*)
1.195 -fun uv_mod_null (p1:int list,0) = p1
1.196 - | uv_mod_null (p1:int list,n1:int) = uv_mod_null(p1,n1-1) @ [0];
1.197 -
1.198 -(*. uv polynomial division, result is (quotient, remainder) .*)
1.199 -(*. only for uv_mod_divides .*)
1.200 -(* FIXME: Division von x^9+x^5+1 durch x-1000 funktioniert nicht integer zu klein *)
1.201 -fun uv_mod_pdiv (p1:uv_poly) ([]:uv_poly) = raise error ("RATIONALS_UV_MOD_PDIV_EXCEPTION: division by zero")
1.202 - | uv_mod_pdiv p1 [x] =
1.203 - let
1.204 - val xs=ref [];
1.205 - in
1.206 - if x<>0 then
1.207 - (
1.208 - xs:=(uv_mod_rem_poly(p1,x));
1.209 - while length(!xs)>0 andalso hd(!xs)=0 do xs:=tl(!xs)
1.210 - )
1.211 - else raise error ("RATIONALS_UV_MOD_PDIV_EXCEPTION: division by zero");
1.212 - ([]:uv_poly,!xs:uv_poly)
1.213 - end
1.214 - | uv_mod_pdiv p1 p2 =
1.215 - let
1.216 - val n= uv_mod_deg(p2);
1.217 - val m= ref (uv_mod_deg(p1));
1.218 - val p1'=ref (rev(p1));
1.219 - val p2'=(rev(p2));
1.220 - val lc2=hd(p2');
1.221 - val q=ref [];
1.222 - val c=ref 0;
1.223 - val output=ref ([],[]);
1.224 - in
1.225 - (
1.226 - if (!m)=0 orelse p2=[0] then raise error ("RATIONALS_UV_MOD_PDIV_EXCEPTION: Division by zero")
1.227 - else
1.228 - (
1.229 - if (!m)<n then
1.230 - (
1.231 - output:=([0],p1)
1.232 - )
1.233 - else
1.234 - (
1.235 - while (!m)>=n do
1.236 - (
1.237 - c:=hd(!p1') div hd(p2');
1.238 - if !c<>0 then
1.239 - (
1.240 - p1':=uv_mod_add_poly(!p1',uv_mod_null(uv_mod_smul_poly(p2',~(!c)),!m-n));
1.241 - while length(!p1')>0 andalso hd(!p1')=0 do p1':= tl(!p1');
1.242 - m:=uv_mod_deg(!p1')
1.243 - )
1.244 - else m:=0
1.245 - );
1.246 - output:=(rev(!q),rev(!p1'))
1.247 - )
1.248 - );
1.249 - !output
1.250 - )
1.251 - end;
1.252 -
1.253 -(*. divides p1 by p2 in Zp .*)
1.254 -fun uv_mod_pdivp (p1:uv_poly) (p2:uv_poly) p =
1.255 - let
1.256 - val n=uv_mod_deg(p2);
1.257 - val m=ref (uv_mod_deg(uv_mod_list_modp p1 p));
1.258 - val p1'=ref (rev(p1));
1.259 - val p2'=(rev(uv_mod_list_modp p2 p));
1.260 - val lc2=hd(p2');
1.261 - val q=ref [];
1.262 - val c=ref 0;
1.263 - val output=ref ([],[]);
1.264 - in
1.265 - (
1.266 - if (!m)=0 orelse p2=[0] then raise error ("RATIONALS_UV_MOD_PDIVP_EXCEPTION: Division by zero")
1.267 - else
1.268 - (
1.269 - if (!m)<n then
1.270 - (
1.271 - output:=([0],p1)
1.272 - )
1.273 - else
1.274 - (
1.275 - while (!m)>=n do
1.276 - (
1.277 - c:=uv_mod_mod2(hd(!p1')*(power lc2 1), p);
1.278 - q:=(!c)::(!q);
1.279 - p1':=uv_mod_list_modp(tl(uv_mod_add_poly(uv_mod_smul_poly(!p1',lc2),
1.280 - uv_mod_smul_poly(uv_mod_smul_poly(p2',hd(!p1')),~1)))) p;
1.281 - m:=(!m)-1
1.282 - );
1.283 -
1.284 - while !p1'<>[] andalso hd(!p1')=0 do
1.285 - (
1.286 - p1':=tl(!p1')
1.287 - );
1.288 -
1.289 - output:=(rev(uv_mod_list_modp (!q) (p)),rev(!p1'))
1.290 - )
1.291 - );
1.292 - !output:uv_poly * uv_poly
1.293 - )
1.294 - end;
1.295 -
1.296 -(*. calculates the remainder of p1/p2 .*)
1.297 -fun uv_mod_prest (p1:uv_poly) ([]:uv_poly) = raise error("UV_MOD_PREST_EXCEPTION: Division by zero")
1.298 - | uv_mod_prest [] p2 = []:uv_poly
1.299 - | uv_mod_prest p1 p2 = (#2(uv_mod_pdiv p1 p2));
1.300 -
1.301 -(*. calculates the remainder of p1/p2 in Zp .*)
1.302 -fun uv_mod_prestp (p1:uv_poly) ([]:uv_poly) p= raise error("UV_MOD_PRESTP_EXCEPTION: Division by zero")
1.303 - | uv_mod_prestp [] p2 p= []:uv_poly
1.304 - | uv_mod_prestp p1 p2 p = #2(uv_mod_pdivp p1 p2 p);
1.305 -
1.306 -(*. calculates the content of a uv polynomial .*)
1.307 -fun uv_mod_cont ([]:uv_poly) = 0
1.308 - | uv_mod_cont (x::p)= gcd_int x (uv_mod_cont(p));
1.309 -
1.310 -(*. divides each coefficient of a uv polynomial by y .*)
1.311 -fun uv_mod_div_list (p:uv_poly,0) = raise error("UV_MOD_DIV_LIST_EXCEPTION: Division by zero")
1.312 - | uv_mod_div_list ([],y) = []:uv_poly
1.313 - | uv_mod_div_list (x::p,y) = (x div y)::uv_mod_div_list(p,y);
1.314 -
1.315 -(*. calculates the primitiv part of a uv polynomial .*)
1.316 -fun uv_mod_pp ([]:uv_poly) = []:uv_poly
1.317 - | uv_mod_pp p =
1.318 - let
1.319 - val c=ref 0;
1.320 - in
1.321 - (
1.322 - c:=uv_mod_cont(p);
1.323 -
1.324 - if !c=0 then raise error ("RATIONALS_UV_MOD_PP_EXCEPTION: content is 0")
1.325 - else uv_mod_div_list(p,!c)
1.326 - )
1.327 - end;
1.328 -
1.329 -(*. gets the leading coefficient of a uv polynomial .*)
1.330 -fun uv_mod_lc ([]:uv_poly) = 0
1.331 - | uv_mod_lc p = hd(rev(p));
1.332 -
1.333 -(*. calculates the euklidean polynomial remainder sequence in Zp .*)
1.334 -fun uv_mod_prs_euklid_p(p1:uv_poly,p2:uv_poly,p)=
1.335 - let
1.336 - val f =ref [];
1.337 - val f'=ref p2;
1.338 - val fi=ref [];
1.339 - in
1.340 - (
1.341 - f:=p2::p1::[];
1.342 - while uv_mod_deg(!f')>0 do
1.343 - (
1.344 - f':=uv_mod_prestp (hd(tl(!f))) (hd(!f)) p;
1.345 - if (!f')<>[] then
1.346 - (
1.347 - fi:=(!f');
1.348 - f:=(!fi)::(!f)
1.349 - )
1.350 - else ()
1.351 - );
1.352 - (!f)
1.353 -
1.354 - )
1.355 - end;
1.356 -
1.357 -(*. calculates the gcd of p1 and p2 in Zp .*)
1.358 -fun uv_mod_gcd_modp ([]:uv_poly) (p2:uv_poly) p = p2:uv_poly
1.359 - | uv_mod_gcd_modp p1 [] p= p1
1.360 - | uv_mod_gcd_modp p1 p2 p=
1.361 - let
1.362 - val p1'=ref[];
1.363 - val p2'=ref[];
1.364 - val pc=ref[];
1.365 - val g=ref [];
1.366 - val d=ref 0;
1.367 - val prs=ref [];
1.368 - in
1.369 - (
1.370 - if uv_mod_deg(p1)>=uv_mod_deg(p2) then
1.371 - (
1.372 - p1':=uv_mod_list_modp (uv_mod_pp(p1)) p;
1.373 - p2':=uv_mod_list_modp (uv_mod_pp(p2)) p
1.374 - )
1.375 - else
1.376 - (
1.377 - p1':=uv_mod_list_modp (uv_mod_pp(p2)) p;
1.378 - p2':=uv_mod_list_modp (uv_mod_pp(p1)) p
1.379 - );
1.380 - d:=uv_mod_mod2((gcd_int (uv_mod_cont(p1))) (uv_mod_cont(p2)), p) ;
1.381 - if !d>(p div 2) then d:=(!d)-p else ();
1.382 -
1.383 - prs:=uv_mod_prs_euklid_p(!p1',!p2',p);
1.384 -
1.385 - if hd(!prs)=[] then pc:=hd(tl(!prs))
1.386 - else pc:=hd(!prs);
1.387 -
1.388 - g:=uv_mod_smul_poly(uv_mod_pp(!pc),!d);
1.389 - !g
1.390 - )
1.391 - end;
1.392 -
1.393 -(*. calculates the minimum of two real values x and y .*)
1.394 -fun uv_mod_r_min(x,y):BasisLibrary.Real.real = if x>y then y else x;
1.395 -
1.396 -(*. calculates the minimum of two integer values x and y .*)
1.397 -fun uv_mod_min(x,y) = if x>y then y else x;
1.398 -
1.399 -(*. adds the squared values of a integer list .*)
1.400 -fun uv_mod_add_qu [] = 0.0
1.401 - | uv_mod_add_qu (x::p) = BasisLibrary.Real.fromInt(x)*BasisLibrary.Real.fromInt(x) + uv_mod_add_qu p;
1.402 -
1.403 -(*. calculates the euklidean norm .*)
1.404 -fun uv_mod_norm ([]:uv_poly) = 0.0
1.405 - | uv_mod_norm p = Math.sqrt(uv_mod_add_qu(p));
1.406 -
1.407 -(*. multipies two values a and b .*)
1.408 -fun uv_mod_multi a b = a * b;
1.409 -
1.410 -(*. decides if x is a prim, the list contains all primes which are lower then x .*)
1.411 -fun uv_mod_prim(x,[])= false
1.412 - | uv_mod_prim(x,[y])=if ((x mod y) <> 0) then true
1.413 - else false
1.414 - | uv_mod_prim(x,y::ys) = if uv_mod_prim(x,[y])
1.415 - then
1.416 - if uv_mod_prim(x,ys) then true
1.417 - else false
1.418 - else false;
1.419 -
1.420 -(*. gets the first prime, which is greater than p and does not divide g .*)
1.421 -fun uv_mod_nextprime(g,p)=
1.422 - let
1.423 - val list=ref [2];
1.424 - val exit=ref 0;
1.425 - val i = ref 2
1.426 - in
1.427 - while (!i<p) do (* calculates the primes lower then p *)
1.428 - (
1.429 - if uv_mod_prim(!i,!list) then
1.430 - (
1.431 - if (p mod !i <> 0)
1.432 - then
1.433 - (
1.434 - list:= (!i)::(!list);
1.435 - i:= (!i)+1
1.436 - )
1.437 - else i:=(!i)+1
1.438 - )
1.439 - else i:= (!i)+1
1.440 - );
1.441 - i:=(p+1);
1.442 - while (!exit=0) do (* calculate next prime which does not divide g *)
1.443 - (
1.444 - if uv_mod_prim(!i,!list) then
1.445 - (
1.446 - if (g mod !i <> 0)
1.447 - then
1.448 - (
1.449 - list:= (!i)::(!list);
1.450 - exit:= (!i)
1.451 - )
1.452 - else i:=(!i)+1
1.453 - )
1.454 - else i:= (!i)+1
1.455 - );
1.456 - !exit
1.457 - end;
1.458 -
1.459 -(*. decides if p1 is a factor of p2 in Zp .*)
1.460 -fun uv_mod_dividesp ([]:uv_poly) (p2:uv_poly) p= raise error("UV_MOD_DIVIDESP: Division by zero")
1.461 - | uv_mod_dividesp p1 p2 p= if uv_mod_prestp p2 p1 p = [] then true else false;
1.462 -
1.463 -(*. decides if p1 is a factor of p2 .*)
1.464 -fun uv_mod_divides ([]:uv_poly) (p2:uv_poly) = raise error("UV_MOD_DIVIDES: Division by zero")
1.465 - | uv_mod_divides p1 p2 = if uv_mod_prest p2 p1 = [] then true else false;
1.466 -
1.467 -(*. chinese remainder algorithm .*)
1.468 -fun uv_mod_cra2(r1,r2,m1,m2)=
1.469 - let
1.470 - val c=ref 0;
1.471 - val r1'=ref 0;
1.472 - val d=ref 0;
1.473 - val a=ref 0;
1.474 - in
1.475 - (
1.476 - while (uv_mod_mod2((!c)*m1,m2))<>1 do
1.477 - (
1.478 - c:=(!c)+1
1.479 - );
1.480 - r1':= uv_mod_mod2(r1,m1);
1.481 - d:=uv_mod_mod2(((r2-(!r1'))*(!c)),m2);
1.482 - !r1'+(!d)*m1
1.483 - )
1.484 - end;
1.485 -
1.486 -(*. applies the chinese remainder algorithmen to the coefficients of x1 and x2 .*)
1.487 -fun uv_mod_cra_2 ([],[],m1,m2) = []
1.488 - | uv_mod_cra_2 ([],x2,m1,m2) = raise error("UV_MOD_CRA_2_EXCEPTION: invalid call x1")
1.489 - | uv_mod_cra_2 (x1,[],m1,m2) = raise error("UV_MOD_CRA_2_EXCEPTION: invalid call x2")
1.490 - | 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));
1.491 -
1.492 -(*. calculates the gcd of two uv polynomials p1' and p2' with the modular algorithm .*)
1.493 -fun uv_mod_gcd (p1':uv_poly) (p2':uv_poly) =
1.494 - let
1.495 - val p1=ref (uv_mod_pp(p1'));
1.496 - val p2=ref (uv_mod_pp(p2'));
1.497 - val c=gcd_int (uv_mod_cont(p1')) (uv_mod_cont(p2'));
1.498 - val temp=ref [];
1.499 - val cp=ref [];
1.500 - val qp=ref [];
1.501 - val q=ref[];
1.502 - val pn=ref 0;
1.503 - val d=ref 0;
1.504 - val g1=ref 0;
1.505 - val p=ref 0;
1.506 - val m=ref 0;
1.507 - val exit=ref 0;
1.508 - val i=ref 1;
1.509 - in
1.510 - if length(!p1)>length(!p2) then ()
1.511 - else
1.512 - (
1.513 - temp:= !p1;
1.514 - p1:= !p2;
1.515 - p2:= !temp
1.516 - );
1.517 -
1.518 -
1.519 - d:=gcd_int (uv_mod_lc(!p1)) (uv_mod_lc(!p2));
1.520 - g1:=uv_mod_lc(!p1)*uv_mod_lc(!p2);
1.521 - p:=4;
1.522 -
1.523 - m:=BasisLibrary.Real.ceil(2.0 *
1.524 - BasisLibrary.Real.fromInt(!d) *
1.525 - BasisLibrary.Real.fromInt(power 2 (uv_mod_min(uv_mod_deg(!p2),uv_mod_deg(!p1)))) *
1.526 - BasisLibrary.Real.fromInt(!d) *
1.527 - uv_mod_r_min(uv_mod_norm(!p1) / BasisLibrary.Real.fromInt(abs(uv_mod_lc(!p1))),
1.528 - uv_mod_norm(!p2) / BasisLibrary.Real.fromInt(abs(uv_mod_lc(!p2)))));
1.529 -
1.530 - while (!exit=0) do
1.531 - (
1.532 - p:=uv_mod_nextprime(!d,!p);
1.533 - cp:=(uv_mod_gcd_modp (uv_mod_list_modp(!p1) (!p)) (uv_mod_list_modp(!p2) (!p)) (!p)) ;
1.534 - if abs(uv_mod_lc(!cp))<>1 then (* leading coefficient = 1 ? *)
1.535 - (
1.536 - i:=1;
1.537 - while (!i)<(!p) andalso (abs(uv_mod_mod2((uv_mod_lc(!cp)*(!i)),(!p)))<>1) do
1.538 - (
1.539 - i:=(!i)+1
1.540 - );
1.541 - cp:=uv_mod_list_modp (map (uv_mod_multi (!i)) (!cp)) (!p)
1.542 - )
1.543 - else ();
1.544 -
1.545 - qp:= ((map (uv_mod_multi (uv_mod_mod2(!d,!p)))) (!cp));
1.546 -
1.547 - if uv_mod_deg(!qp)=0 then (q:=[1]; exit:=1) else ();
1.548 -
1.549 - pn:=(!p);
1.550 - q:=(!qp);
1.551 -
1.552 - while !pn<= !m andalso !m>(!p) andalso !exit=0 do
1.553 - (
1.554 - p:=uv_mod_nextprime(!d,!p);
1.555 - cp:=(uv_mod_gcd_modp (uv_mod_list_modp(!p1) (!p)) (uv_mod_list_modp(!p2) (!p)) (!p));
1.556 - if uv_mod_lc(!cp)<>1 then (* leading coefficient = 1 ? *)
1.557 - (
1.558 - i:=1;
1.559 - while (!i)<(!p) andalso ((uv_mod_mod2((uv_mod_lc(!q)*(!i)),(!p)))<>1) do
1.560 - (
1.561 - i:=(!i)+1
1.562 - );
1.563 - cp:=uv_mod_list_modp (map (uv_mod_multi (!i)) (!cp)) (!p)
1.564 - )
1.565 - else ();
1.566 -
1.567 - qp:=uv_mod_list_modp ((map (uv_mod_multi (uv_mod_mod2(!d,!p)))) (!cp) ) (!p);
1.568 - if uv_mod_deg(!qp)>uv_mod_deg(!q) then
1.569 - (
1.570 - (*print("degree to high!!!\n")*)
1.571 - )
1.572 - else
1.573 - (
1.574 - if uv_mod_deg(!qp)=uv_mod_deg(!q) then
1.575 - (
1.576 - q:=uv_mod_cra_2(!q,!qp,!pn,!p);
1.577 - pn:=(!pn) * !p;
1.578 - q:=uv_mod_pp(uv_mod_list_modp (!q) (!pn)); (* found already gcd ? *)
1.579 - if (uv_mod_divides (!q) (p1')) andalso (uv_mod_divides (!q) (p2')) then (exit:=1) else ()
1.580 - )
1.581 - else
1.582 - (
1.583 - if uv_mod_deg(!qp)<uv_mod_deg(!q) then
1.584 - (
1.585 - pn:= !p;
1.586 - q:= !qp
1.587 - )
1.588 - else ()
1.589 - )
1.590 - )
1.591 - );
1.592 - q:=uv_mod_pp(uv_mod_list_modp (!q) (!pn));
1.593 - if (uv_mod_divides (!q) (p1')) andalso (uv_mod_divides (!q) (p2')) then exit:=1 else ()
1.594 - );
1.595 - uv_mod_smul_poly(!q,c):uv_poly
1.596 - end;
1.597 -
1.598 -(*. multivariate polynomials .*)
1.599 -(*. multivariate polynomials are represented as a list of the pairs,
1.600 - first is the coefficent and the second is a list of the exponents .*)
1.601 -(*. 5 * x^5 * y^3 + 4 * x^3 * z^2 + 2 * x^2 * y * z^3 - z - 19
1.602 - => [(5,[5,3,0]),(4,[3,0,2]),(2,[2,1,3]),(~1,[0,0,1]),(~19,[0,0,0])] .*)
1.603 -
1.604 -(*. global variables .*)
1.605 -(*. order indicators .*)
1.606 -val LEX_=0; (* lexicographical term order *)
1.607 -val GGO_=1; (* greatest degree order *)
1.608 -
1.609 -(*. datatypes for internal representation.*)
1.610 -type mv_monom = (int * (*.coefficient or the monom.*)
1.611 - int list); (*.list of exponents) .*)
1.612 -fun mv_monom2str (i, is) = "("^ int2str i^"," ^ ints2str' is ^ ")";
1.613 -
1.614 -type mv_poly = mv_monom list;
1.615 -fun mv_poly2str p = (strs2str' o (map mv_monom2str)) p;
1.616 -
1.617 -(*. help function for monom_greater and geq .*)
1.618 -fun mv_mg_hlp([]) = EQUAL
1.619 - | mv_mg_hlp(x::list)=if x<0 then LESS
1.620 - else if x>0 then GREATER
1.621 - else mv_mg_hlp(list);
1.622 -
1.623 -(*. adds a list of values .*)
1.624 -fun mv_addlist([]) = 0
1.625 - | mv_addlist(p1) = hd(p1)+mv_addlist(tl(p1));
1.626 -
1.627 -(*. tests if the monomial M1 is greater as the monomial M2 and returns a boolean value .*)
1.628 -(*. 2 orders are implemented LEX_/GGO_ (lexigraphical/greatest degree order) .*)
1.629 -fun mv_monom_greater((M1x,M1l):mv_monom,(M2x,M2l):mv_monom,order)=
1.630 - if order=LEX_ then
1.631 - (
1.632 - if length(M1l)<>length(M2l) then raise error ("RATIONALS_MV_MONOM_GREATER_EXCEPTION: Order error")
1.633 - else if (mv_mg_hlp((map op- (M1l~~M2l)))<>GREATER) then false else true
1.634 - )
1.635 - else
1.636 - if order=GGO_ then
1.637 - (
1.638 - if length(M1l)<>length(M2l) then raise error ("RATIONALS_MV_MONOM_GREATER_EXCEPTION: Order error")
1.639 - else
1.640 - if mv_addlist(M1l)=mv_addlist(M2l) then if (mv_mg_hlp((map op- (M1l~~M2l)))<>GREATER) then false else true
1.641 - else if mv_addlist(M1l)>mv_addlist(M2l) then true else false
1.642 - )
1.643 - else raise error ("RATIONALS_MV_MONOM_GREATER_EXCEPTION: Wrong Order");
1.644 -
1.645 -(*. tests if the monomial X is greater as the monomial Y and returns a order value (GREATER,EQUAL,LESS) .*)
1.646 -(*. 2 orders are implemented LEX_/GGO_ (lexigraphical/greatest degree order) .*)
1.647 -fun mv_geq order ((x1,x):mv_monom,(x2,y):mv_monom) =
1.648 -let
1.649 - val temp=ref EQUAL;
1.650 -in
1.651 - if order=LEX_ then
1.652 - (
1.653 - if length(x)<>length(y) then
1.654 - raise error ("RATIONALS_MV_GEQ_EXCEPTION: Order error")
1.655 - else
1.656 - (
1.657 - temp:=mv_mg_hlp((map op- (x~~y)));
1.658 - if !temp=EQUAL then
1.659 - ( if x1=x2 then EQUAL
1.660 - else if x1>x2 then GREATER
1.661 - else LESS
1.662 - )
1.663 - else (!temp)
1.664 - )
1.665 - )
1.666 - else
1.667 - if order=GGO_ then
1.668 - (
1.669 - if length(x)<>length(y) then
1.670 - raise error ("RATIONALS_MV_GEQ_EXCEPTION: Order error")
1.671 - else
1.672 - if mv_addlist(x)=mv_addlist(y) then
1.673 - (mv_mg_hlp((map op- (x~~y))))
1.674 - else if mv_addlist(x)>mv_addlist(y) then GREATER else LESS
1.675 - )
1.676 - else raise error ("RATIONALS_MV_GEQ_EXCEPTION: Wrong Order")
1.677 -end;
1.678 -
1.679 -(*. cuts the first variable from a polynomial .*)
1.680 -fun mv_cut([]:mv_poly)=[]:mv_poly
1.681 - | mv_cut((x,[])::list) = raise error ("RATIONALS_MV_CUT_EXCEPTION: Invalid list ")
1.682 - | mv_cut((x,y::ys)::list)=(x,ys)::mv_cut(list);
1.683 -
1.684 -(*. leading power product .*)
1.685 -fun mv_lpp([]:mv_poly,order) = []
1.686 - | mv_lpp([(x,y)],order) = y
1.687 - | mv_lpp(p1,order) = #2(hd(rev(sort (mv_geq order) p1)));
1.688 -
1.689 -(*. leading monomial .*)
1.690 -fun mv_lm([]:mv_poly,order) = (0,[]):mv_monom
1.691 - | mv_lm([x],order) = x
1.692 - | mv_lm(p1,order) = hd(rev(sort (mv_geq order) p1));
1.693 -
1.694 -(*. leading coefficient in term order .*)
1.695 -fun mv_lc2([]:mv_poly,order) = 0
1.696 - | mv_lc2([(x,y)],order) = x
1.697 - | mv_lc2(p1,order) = #1(hd(rev(sort (mv_geq order) p1)));
1.698 -
1.699 -
1.700 -(*. reverse the coefficients in mv polynomial .*)
1.701 -fun mv_rev_to([]:mv_poly) = []:mv_poly
1.702 - | mv_rev_to((c,e)::xs) = (c,rev(e))::mv_rev_to(xs);
1.703 -
1.704 -(*. leading coefficient in reverse term order .*)
1.705 -fun mv_lc([]:mv_poly,order) = []:mv_poly
1.706 - | mv_lc([(x,y)],order) = mv_rev_to(mv_cut(mv_rev_to([(x,y)])))
1.707 - | mv_lc(p1,order) =
1.708 - let
1.709 - val p1o=ref (rev(sort (mv_geq order) (mv_rev_to(p1))));
1.710 - val lp=hd(#2(hd(!p1o)));
1.711 - val lc=ref [];
1.712 - in
1.713 - (
1.714 - while (length(!p1o)>0 andalso hd(#2(hd(!p1o)))=lp) do
1.715 - (
1.716 - lc:=hd(mv_cut([hd(!p1o)]))::(!lc);
1.717 - p1o:=tl(!p1o)
1.718 - );
1.719 - if !lc=[] then raise error ("RATIONALS_MV_LC_EXCEPTION: lc is empty") else ();
1.720 - mv_rev_to(!lc)
1.721 - )
1.722 - end;
1.723 -
1.724 -(*. compares two powerproducts .*)
1.725 -fun mv_monom_equal((_,xlist):mv_monom,(_,ylist):mv_monom) = (foldr and_) (((map op=) (xlist~~ylist)),true);
1.726 -
1.727 -(*. help function for mv_add .*)
1.728 -fun mv_madd([]:mv_poly,[]:mv_poly,order) = []:mv_poly
1.729 - | mv_madd([(0,_)],p2,order) = p2
1.730 - | mv_madd(p1,[(0,_)],order) = p1
1.731 - | mv_madd([],p2,order) = p2
1.732 - | mv_madd(p1,[],order) = p1
1.733 - | mv_madd(p1,p2,order) =
1.734 - (
1.735 - if mv_monom_greater(hd(p1),hd(p2),order)
1.736 - then hd(p1)::mv_madd(tl(p1),p2,order)
1.737 - else if mv_monom_equal(hd(p1),hd(p2))
1.738 - then if mv_lc2(p1,order)+mv_lc2(p2,order)<>0
1.739 - then (mv_lc2(p1,order)+mv_lc2(p2,order),mv_lpp(p1,order))::mv_madd(tl(p1),tl(p2),order)
1.740 - else mv_madd(tl(p1),tl(p2),order)
1.741 - else hd(p2)::mv_madd(p1,tl(p2),order)
1.742 - )
1.743 -
1.744 -(*. adds two multivariate polynomials .*)
1.745 -fun mv_add([]:mv_poly,p2:mv_poly,order) = p2
1.746 - | mv_add(p1,[],order) = p1
1.747 - | mv_add(p1,p2,order) = mv_madd(rev(sort (mv_geq order) p1),rev(sort (mv_geq order) p2), order);
1.748 -
1.749 -(*. monom multiplication .*)
1.750 -fun mv_mmul((x1,y1):mv_monom,(x2,y2):mv_monom)=(x1*x2,(map op+) (y1~~y2)):mv_monom;
1.751 -
1.752 -(*. deletes all monomials with coefficient 0 .*)
1.753 -fun mv_shorten([]:mv_poly,order) = []:mv_poly
1.754 - | mv_shorten(x::xs,order)=mv_madd([x],mv_shorten(xs,order),order);
1.755 -
1.756 -(*. zeros a list .*)
1.757 -fun mv_null2([])=[]
1.758 - | mv_null2(x::l)=0::mv_null2(l);
1.759 -
1.760 -(*. multiplies two multivariate polynomials .*)
1.761 -fun mv_mul([]:mv_poly,[]:mv_poly,_) = []:mv_poly
1.762 - | mv_mul([],y::p2,_) = [(0,mv_null2(#2(y)))]
1.763 - | mv_mul(x::p1,[],_) = [(0,mv_null2(#2(x)))]
1.764 - | mv_mul(x::p1,y::p2,order) = mv_shorten(rev(sort (mv_geq order) (mv_mmul(x,y) :: (mv_mul(p1,y::p2,order) @
1.765 - mv_mul([x],p2,order)))),order);
1.766 -
1.767 -(*. gets the maximum value of a list .*)
1.768 -fun mv_getmax([])=0
1.769 - | mv_getmax(x::p1)= let
1.770 - val m=mv_getmax(p1);
1.771 - in
1.772 - if m>x then m
1.773 - else x
1.774 - end;
1.775 -(*. calculates the maximum degree of an multivariate polynomial .*)
1.776 -fun mv_grad([]:mv_poly) = 0
1.777 - | mv_grad(p1:mv_poly)= mv_getmax((map mv_addlist) ((map #2) p1));
1.778 -
1.779 -(*. converts the sign of a value .*)
1.780 -fun mv_minus(x)=(~1) * x;
1.781 -
1.782 -(*. converts the sign of all coefficients of a polynomial .*)
1.783 -fun mv_minus2([]:mv_poly)=[]:mv_poly
1.784 - | mv_minus2(p1)=(mv_minus(#1(hd(p1))),#2(hd(p1)))::(mv_minus2(tl(p1)));
1.785 -
1.786 -(*. searches for a negativ value in a list .*)
1.787 -fun mv_is_negativ([])=false
1.788 - | mv_is_negativ(x::xs)=if x<0 then true else mv_is_negativ(xs);
1.789 -
1.790 -(*. division of monomials .*)
1.791 -fun mv_mdiv((0,[]):mv_monom,_:mv_monom)=(0,[]):mv_monom
1.792 - | mv_mdiv(_,(0,[]))= raise error ("RATIONALS_MV_MDIV_EXCEPTION Division by 0 ")
1.793 - | mv_mdiv(p1:mv_monom,p2:mv_monom)=
1.794 - let
1.795 - val c=ref (#1(p2));
1.796 - val pp=ref [];
1.797 - in
1.798 - (
1.799 - if !c=0 then raise error("MV_MDIV_EXCEPTION Dividing by zero")
1.800 - else c:=(#1(p1) div #1(p2));
1.801 - if #1(p2)<>0 then
1.802 - (
1.803 - pp:=(#2(mv_mmul((1,#2(p1)),(1,(map mv_minus) (#2(p2))))));
1.804 - if mv_is_negativ(!pp) then (0,!pp)
1.805 - else (!c,!pp)
1.806 - )
1.807 - else raise error("MV_MDIV_EXCEPTION Dividing by empty Polynom")
1.808 - )
1.809 - end;
1.810 -
1.811 -(*. prints a polynom for (internal use only) .*)
1.812 -fun mv_print_poly([]:mv_poly)=print("[]\n")
1.813 - | mv_print_poly((x,y)::[])= print("("^BasisLibrary.Int.toString(x)^","^ints2str(y)^")\n")
1.814 - | mv_print_poly((x,y)::p1) = (print("("^BasisLibrary.Int.toString(x)^","^ints2str(y)^"),");mv_print_poly(p1));
1.815 -
1.816 -
1.817 -(*. division of two multivariate polynomials .*)
1.818 -fun mv_division([]:mv_poly,g:mv_poly,order)=([]:mv_poly,[]:mv_poly)
1.819 - | mv_division(f,[],order)= raise error ("RATIONALS_MV_DIVISION_EXCEPTION Division by zero")
1.820 - | mv_division(f,g,order)=
1.821 - let
1.822 - val r=ref [];
1.823 - val q=ref [];
1.824 - val g'=ref [];
1.825 - val k=ref 0;
1.826 - val m=ref (0,[0]);
1.827 - val exit=ref 0;
1.828 - in
1.829 - r := rev(sort (mv_geq order) (mv_shorten(f,order)));
1.830 - g':= rev(sort (mv_geq order) (mv_shorten(g,order)));
1.831 - if #1(hd(!g'))=0 then raise error("RATIONALS_MV_DIVISION_EXCEPTION: Dividing by zero") else ();
1.832 - if (mv_monom_greater (hd(!g'),hd(!r),order)) then ([(0,mv_null2(#2(hd(f))))],(!r))
1.833 - else
1.834 - (
1.835 - exit:=0;
1.836 - while (if (!exit)=0 then not(mv_monom_greater (hd(!g'),hd(!r),order)) else false) do
1.837 - (
1.838 - if (#1(mv_lm(!g',order)))<>0 then m:=mv_mdiv(mv_lm(!r,order),mv_lm(!g',order))
1.839 - else raise error ("RATIONALS_MV_DIVISION_EXCEPTION: Dividing by zero");
1.840 - if #1(!m)<>0 then
1.841 - (
1.842 - q:=(!m)::(!q);
1.843 - r:=mv_add((!r),mv_minus2(mv_mul(!g',[!m],order)),order)
1.844 - )
1.845 - else exit:=1;
1.846 - if (if length(!r)<>0 then length(!g')<>0 else false) then ()
1.847 - else (exit:=1)
1.848 - );
1.849 - (rev(!q),!r)
1.850 - )
1.851 - end;
1.852 -
1.853 -(*. multiplies a polynomial with an integer .*)
1.854 -fun mv_skalar_mul([]:mv_poly,c) = []:mv_poly
1.855 - | mv_skalar_mul((x,y)::p1,c) = ((x * c),y)::mv_skalar_mul(p1,c);
1.856 -
1.857 -(*. inserts the a first variable into an polynomial with exponent v .*)
1.858 -fun mv_correct([]:mv_poly,v:int)=[]:mv_poly
1.859 - | mv_correct((x,y)::list,v:int)=(x,v::y)::mv_correct(list,v);
1.860 -
1.861 -(*. multivariate case .*)
1.862 -
1.863 -(*. decides if x is a factor of y .*)
1.864 -fun mv_divides([]:mv_poly,[]:mv_poly)= raise error("RATIONALS_MV_DIVIDES_EXCEPTION: division by zero")
1.865 - | mv_divides(x,[]) = raise error("RATIONALS_MV_DIVIDES_EXCEPTION: division by zero")
1.866 - | mv_divides(x:mv_poly,y:mv_poly) = #2(mv_division(y,x,LEX_))=[];
1.867 -
1.868 -(*. gets the maximum of a and b .*)
1.869 -fun mv_max(a,b) = if a>b then a else b;
1.870 -
1.871 -(*. gets the maximum exponent of a mv polynomial in the lexicographic term order .*)
1.872 -fun mv_deg([]:mv_poly) = 0
1.873 - | mv_deg(p1)=
1.874 - let
1.875 - val p1'=mv_shorten(p1,LEX_);
1.876 - in
1.877 - if length(p1')=0 then 0
1.878 - else mv_max(hd(#2(hd(p1'))),mv_deg(tl(p1')))
1.879 - end;
1.880 -
1.881 -(*. gets the maximum exponent of a mv polynomial in the reverse lexicographic term order .*)
1.882 -fun mv_deg2([]:mv_poly) = 0
1.883 - | mv_deg2(p1)=
1.884 - let
1.885 - val p1'=mv_shorten(p1,LEX_);
1.886 - in
1.887 - if length(p1')=0 then 0
1.888 - else mv_max(hd(rev(#2(hd(p1')))),mv_deg2(tl(p1')))
1.889 - end;
1.890 -
1.891 -(*. evaluates the mv polynomial at the value v of the main variable .*)
1.892 -fun mv_subs([]:mv_poly,v) = []:mv_poly
1.893 - | mv_subs((c,e)::p1:mv_poly,v) = mv_skalar_mul(mv_cut([(c,e)]),power v (hd(e))) @ mv_subs(p1,v);
1.894 -
1.895 -(*. calculates the content of a uv-polynomial in mv-representation .*)
1.896 -fun uv_content2([]:mv_poly) = 0
1.897 - | uv_content2((c,e)::p1) = (gcd_int c (uv_content2(p1)));
1.898 -
1.899 -(*. converts a uv-polynomial from mv-representation to uv-representation .*)
1.900 -fun uv_to_list ([]:mv_poly)=[]:uv_poly
1.901 - | uv_to_list ((c1,e1)::others) =
1.902 - let
1.903 - val count=ref 0;
1.904 - val max=mv_grad((c1,e1)::others);
1.905 - val help=ref ((c1,e1)::others);
1.906 - val list=ref [];
1.907 - in
1.908 - if length(e1)>1 then raise error ("RATIONALS_TO_LIST_EXCEPTION: not univariate")
1.909 - else if length(e1)=0 then [c1]
1.910 - else
1.911 - (
1.912 - count:=0;
1.913 - while (!count)<=max do
1.914 - (
1.915 - if length(!help)>0 andalso hd(#2(hd(!help)))=max-(!count) then
1.916 - (
1.917 - list:=(#1(hd(!help)))::(!list);
1.918 - help:=tl(!help)
1.919 - )
1.920 - else
1.921 - (
1.922 - list:= 0::(!list)
1.923 - );
1.924 - count := (!count) + 1
1.925 - );
1.926 - (!list)
1.927 - )
1.928 - end;
1.929 -
1.930 -(*. converts a uv-polynomial from uv-representation to mv-representation .*)
1.931 -fun uv_to_poly ([]:uv_poly) = []:mv_poly
1.932 - | uv_to_poly p1 =
1.933 - let
1.934 - val count=ref 0;
1.935 - val help=ref p1;
1.936 - val list=ref [];
1.937 - in
1.938 - while length(!help)>0 do
1.939 - (
1.940 - if hd(!help)=0 then ()
1.941 - else list:=(hd(!help),[!count])::(!list);
1.942 - count:=(!count)+1;
1.943 - help:=tl(!help)
1.944 - );
1.945 - (!list)
1.946 - end;
1.947 -
1.948 -(*. univariate gcd calculation from polynomials in multivariate representation .*)
1.949 -fun uv_gcd ([]:mv_poly) p2 = p2
1.950 - | uv_gcd p1 ([]:mv_poly) = p1
1.951 - | uv_gcd p1 [(c,[e])] =
1.952 - let
1.953 - val list=ref (rev(sort (mv_geq LEX_) (mv_shorten(p1,LEX_))));
1.954 - val min=uv_mod_min(e,(hd(#2(hd(rev(!list))))));
1.955 - in
1.956 - [(gcd_int (uv_content2(p1)) c,[min])]
1.957 - end
1.958 - | uv_gcd [(c,[e])] p2 =
1.959 - let
1.960 - val list=ref (rev(sort (mv_geq LEX_) (mv_shorten(p2,LEX_))));
1.961 - val min=uv_mod_min(e,(hd(#2(hd(rev(!list))))));
1.962 - in
1.963 - [(gcd_int (uv_content2(p2)) c,[min])]
1.964 - end
1.965 - | uv_gcd p11 p22 = uv_to_poly(uv_mod_gcd (uv_to_list(mv_shorten(p11,LEX_))) (uv_to_list(mv_shorten(p22,LEX_))));
1.966 -
1.967 -(*. help function for the newton interpolation .*)
1.968 -fun mv_newton_help ([]:mv_poly list,k:int) = []:mv_poly list
1.969 - | mv_newton_help (pl:mv_poly list,k) =
1.970 - let
1.971 - val x=ref (rev(pl));
1.972 - val t=ref [];
1.973 - val y=ref [];
1.974 - val n=ref 1;
1.975 - val n1=ref[];
1.976 - in
1.977 - (
1.978 - while length(!x)>1 do
1.979 - (
1.980 - if length(hd(!x))>0 then n1:=mv_null2(#2(hd(hd(!x))))
1.981 - else if length(hd(tl(!x)))>0 then n1:=mv_null2(#2(hd(hd(tl(!x)))))
1.982 - else n1:=[];
1.983 - t:= #1(mv_division(mv_add(hd(!x),mv_skalar_mul(hd(tl(!x)),~1),LEX_),[(k,!n1)],LEX_));
1.984 - y:=(!t)::(!y);
1.985 - x:=tl(!x)
1.986 - );
1.987 - (!y)
1.988 - )
1.989 - end;
1.990 -
1.991 -(*. help function for the newton interpolation .*)
1.992 -fun mv_newton_add ([]:mv_poly list) t = []:mv_poly
1.993 - | mv_newton_add [x:mv_poly] t = x
1.994 - | mv_newton_add (pl:mv_poly list) t =
1.995 - let
1.996 - val expos=ref [];
1.997 - val pll=ref pl;
1.998 - in
1.999 - (
1.1000 -
1.1001 - while length(!pll)>0 andalso hd(!pll)=[] do
1.1002 - (
1.1003 - pll:=tl(!pll)
1.1004 - );
1.1005 - if length(!pll)>0 then expos:= #2(hd(hd(!pll))) else expos:=[];
1.1006 - mv_add(hd(pl),
1.1007 - mv_mul(
1.1008 - mv_add(mv_correct(mv_cut([(1,mv_null2(!expos))]),1),[(~t,mv_null2(!expos))],LEX_),
1.1009 - mv_newton_add (tl(pl)) (t+1),
1.1010 - LEX_
1.1011 - ),
1.1012 - LEX_)
1.1013 - )
1.1014 - end;
1.1015 -
1.1016 -(*. calculates the newton interpolation with polynomial coefficients .*)
1.1017 -(*. step-depth is 1 and if the result is not an integerpolynomial .*)
1.1018 -(*. this function returns [] .*)
1.1019 -fun mv_newton ([]:(mv_poly) list) = []:mv_poly
1.1020 - | mv_newton ([mp]:(mv_poly) list) = mp:mv_poly
1.1021 - | mv_newton pl =
1.1022 - let
1.1023 - val c=ref pl;
1.1024 - val c1=ref [];
1.1025 - val n=length(pl);
1.1026 - val k=ref 1;
1.1027 - val i=ref n;
1.1028 - val ppl=ref [];
1.1029 - in
1.1030 - c1:=hd(pl)::[];
1.1031 - c:=mv_newton_help(!c,!k);
1.1032 - c1:=(hd(!c))::(!c1);
1.1033 - while(length(!c)>1 andalso !k<n) do
1.1034 - (
1.1035 - k:=(!k)+1;
1.1036 - while length(!c)>0 andalso hd(!c)=[] do c:=tl(!c);
1.1037 - if !c=[] then () else c:=mv_newton_help(!c,!k);
1.1038 - ppl:= !c;
1.1039 - if !c=[] then () else c1:=(hd(!c))::(!c1)
1.1040 - );
1.1041 - while hd(!c1)=[] do c1:=tl(!c1);
1.1042 - c1:=rev(!c1);
1.1043 - ppl:= !c1;
1.1044 - mv_newton_add (!c1) 1
1.1045 - end;
1.1046 -
1.1047 -(*. sets the exponents of the first variable to zero .*)
1.1048 -fun mv_null3([]:mv_poly) = []:mv_poly
1.1049 - | mv_null3((x,y)::xs) = (x,0::tl(y))::mv_null3(xs);
1.1050 -
1.1051 -(*. calculates the minimum exponents of a multivariate polynomial .*)
1.1052 -fun mv_min_pp([]:mv_poly)=[]
1.1053 - | mv_min_pp((c,e)::xs)=
1.1054 - let
1.1055 - val y=ref xs;
1.1056 - val x=ref [];
1.1057 - in
1.1058 - (
1.1059 - x:=e;
1.1060 - while length(!y)>0 do
1.1061 - (
1.1062 - x:=(map uv_mod_min) ((!x) ~~ (#2(hd(!y))));
1.1063 - y:=tl(!y)
1.1064 - );
1.1065 - !x
1.1066 - )
1.1067 - end;
1.1068 -
1.1069 -(*. checks if all elements of the list have value zero .*)
1.1070 -fun list_is_null [] = true
1.1071 - | list_is_null (x::xs) = (x=0 andalso list_is_null(xs));
1.1072 -
1.1073 -(* check if main variable is zero*)
1.1074 -fun main_zero (ms : mv_poly) = (list_is_null o (map (hd o #2))) ms;
1.1075 -
1.1076 -(*. calculates the content of an polynomial .*)
1.1077 -fun mv_content([]:mv_poly) = []:mv_poly
1.1078 - | mv_content(p1) =
1.1079 - let
1.1080 - val list=ref (rev(sort (mv_geq LEX_) (mv_shorten(p1,LEX_))));
1.1081 - val test=ref (hd(#2(hd(!list))));
1.1082 - val result=ref [];
1.1083 - val min=(hd(#2(hd(rev(!list)))));
1.1084 - in
1.1085 - (
1.1086 - if length(!list)>1 then
1.1087 - (
1.1088 - while (if length(!list)>0 then (hd(#2(hd(!list)))=(!test)) else false) do
1.1089 - (
1.1090 - result:=(#1(hd(!list)),tl(#2(hd(!list))))::(!result);
1.1091 -
1.1092 - if length(!list)<1 then list:=[]
1.1093 - else list:=tl(!list)
1.1094 -
1.1095 - );
1.1096 - if length(!list)>0 then
1.1097 - (
1.1098 - list:=mv_gcd (!result) (mv_cut(mv_content(!list)))
1.1099 - )
1.1100 - else list:=(!result);
1.1101 - list:=mv_correct(!list,0);
1.1102 - (!list)
1.1103 - )
1.1104 - else
1.1105 - (
1.1106 - mv_null3(!list)
1.1107 - )
1.1108 - )
1.1109 - end
1.1110 -
1.1111 -(*. calculates the primitiv part of a polynomial .*)
1.1112 -and mv_pp([]:mv_poly) = []:mv_poly
1.1113 - | mv_pp(p1) = let
1.1114 - val cont=ref [];
1.1115 - val pp=ref[];
1.1116 - in
1.1117 - cont:=mv_content(p1);
1.1118 - pp:=(#1(mv_division(p1,!cont,LEX_)));
1.1119 - if !pp=[]
1.1120 - then raise error("RATIONALS_MV_PP_EXCEPTION: Invalid Content ")
1.1121 - else (!pp)
1.1122 - end
1.1123 -
1.1124 -(*. calculates the gcd of two multivariate polynomials with a modular approach .*)
1.1125 -and mv_gcd ([]:mv_poly) ([]:mv_poly) :mv_poly= []:mv_poly
1.1126 - | mv_gcd ([]:mv_poly) (p2) :mv_poly= p2:mv_poly
1.1127 - | mv_gcd (p1:mv_poly) ([]) :mv_poly= p1:mv_poly
1.1128 - | mv_gcd ([(x,xs)]:mv_poly) ([(y,ys)]):mv_poly =
1.1129 - let
1.1130 - val xpoly:mv_poly = [(x,xs)];
1.1131 - val ypoly:mv_poly = [(y,ys)];
1.1132 - in
1.1133 - (
1.1134 - if xs=ys then [((gcd_int x y),xs)]
1.1135 - else [((gcd_int x y),(map uv_mod_min)(xs~~ys))]:mv_poly
1.1136 - )
1.1137 - end
1.1138 - | mv_gcd (p1:mv_poly) ([(y,ys)]) :mv_poly=
1.1139 - (
1.1140 - [(gcd_int (uv_content2(p1)) (y),(map uv_mod_min)(mv_min_pp(p1)~~ys))]:mv_poly
1.1141 - )
1.1142 - | mv_gcd ([(y,ys)]:mv_poly) (p2):mv_poly =
1.1143 - (
1.1144 - [(gcd_int (uv_content2(p2)) (y),(map uv_mod_min)(mv_min_pp(p2)~~ys))]:mv_poly
1.1145 - )
1.1146 - | mv_gcd (p1':mv_poly) (p2':mv_poly):mv_poly=
1.1147 - let
1.1148 - val vc=length(#2(hd(p1')));
1.1149 - val cont =
1.1150 - (
1.1151 - if main_zero(mv_content(p1')) andalso
1.1152 - (main_zero(mv_content(p2'))) then
1.1153 - mv_correct((mv_gcd (mv_cut(mv_content(p1'))) (mv_cut(mv_content(p2')))),0)
1.1154 - else
1.1155 - mv_gcd (mv_content(p1')) (mv_content(p2'))
1.1156 - );
1.1157 - val p1= #1(mv_division(p1',mv_content(p1'),LEX_));
1.1158 - val p2= #1(mv_division(p2',mv_content(p2'),LEX_));
1.1159 - val gcd=ref [];
1.1160 - val candidate=ref [];
1.1161 - val interpolation_list=ref [];
1.1162 - val delta=ref [];
1.1163 - val p1r = ref [];
1.1164 - val p2r = ref [];
1.1165 - val p1r' = ref [];
1.1166 - val p2r' = ref [];
1.1167 - val factor=ref [];
1.1168 - val r=ref 0;
1.1169 - val gcd_r=ref [];
1.1170 - val d=ref 0;
1.1171 - val exit=ref 0;
1.1172 - val current_degree=ref 99999; (*. FIXME: unlimited ! .*)
1.1173 - in
1.1174 - (
1.1175 - if vc<2 then (* areUnivariate(p1',p2') *)
1.1176 - (
1.1177 - gcd:=uv_gcd (mv_shorten(p1',LEX_)) (mv_shorten(p2',LEX_))
1.1178 - )
1.1179 - else
1.1180 - (
1.1181 - while !exit=0 do
1.1182 - (
1.1183 - r:=(!r)+1;
1.1184 - p1r := mv_lc(p1,LEX_);
1.1185 - p2r := mv_lc(p2,LEX_);
1.1186 - if main_zero(!p1r) andalso
1.1187 - main_zero(!p2r)
1.1188 - then
1.1189 - (
1.1190 - delta := mv_correct((mv_gcd (mv_cut (!p1r)) (mv_cut (!p2r))),0)
1.1191 - )
1.1192 - else
1.1193 - (
1.1194 - delta := mv_gcd (!p1r) (!p2r)
1.1195 - );
1.1196 - (*if mv_shorten(mv_subs(!p1r,!r),LEX_)=[] andalso
1.1197 - mv_shorten(mv_subs(!p2r,!r),LEX_)=[] *)
1.1198 - if mv_lc2(mv_shorten(mv_subs(!p1r,!r),LEX_),LEX_)=0 andalso
1.1199 - mv_lc2(mv_shorten(mv_subs(!p2r,!r),LEX_),LEX_)=0
1.1200 - then
1.1201 - (
1.1202 - )
1.1203 - else
1.1204 - (
1.1205 - gcd_r:=mv_shorten(mv_gcd (mv_shorten(mv_subs(p1,!r),LEX_))
1.1206 - (mv_shorten(mv_subs(p2,!r),LEX_)) ,LEX_);
1.1207 - gcd_r:= #1(mv_division(mv_mul(mv_correct(mv_subs(!delta,!r),0),!gcd_r,LEX_),
1.1208 - mv_correct(mv_lc(!gcd_r,LEX_),0),LEX_));
1.1209 - d:=mv_deg2(!gcd_r); (* deg(gcd_r,z) *)
1.1210 - if (!d < !current_degree) then
1.1211 - (
1.1212 - current_degree:= !d;
1.1213 - interpolation_list:=mv_correct(!gcd_r,0)::(!interpolation_list)
1.1214 - )
1.1215 - else
1.1216 - (
1.1217 - if (!d = !current_degree) then
1.1218 - (
1.1219 - interpolation_list:=mv_correct(!gcd_r,0)::(!interpolation_list)
1.1220 - )
1.1221 - else ()
1.1222 - )
1.1223 - );
1.1224 - if length(!interpolation_list)> uv_mod_min(mv_deg(p1),mv_deg(p2)) then
1.1225 - (
1.1226 - candidate := mv_newton(rev(!interpolation_list));
1.1227 - if !candidate=[] then ()
1.1228 - else
1.1229 - (
1.1230 - candidate:=mv_pp(!candidate);
1.1231 - if mv_divides(!candidate,p1) andalso mv_divides(!candidate,p2) then
1.1232 - (
1.1233 - gcd:= mv_mul(!candidate,cont,LEX_);
1.1234 - exit:=1
1.1235 - )
1.1236 - else ()
1.1237 - );
1.1238 - interpolation_list:=[mv_correct(!gcd_r,0)]
1.1239 - )
1.1240 - else ()
1.1241 - )
1.1242 - );
1.1243 - (!gcd):mv_poly
1.1244 - )
1.1245 - end;
1.1246 -
1.1247 -
1.1248 -(*. calculates the least common divisor of two polynomials .*)
1.1249 -fun mv_lcm (p1:mv_poly) (p2:mv_poly) :mv_poly =
1.1250 - (
1.1251 - #1(mv_division(mv_mul(p1,p2,LEX_),mv_gcd p1 p2,LEX_))
1.1252 - );
1.1253 -
1.1254 -(*. gets the variables (strings) of a term .*)
1.1255 -fun get_vars(term1) = (map free2str) (vars term1); (*["a","b","c"]; *)
1.1256 -
1.1257 -(*. counts the negative coefficents in a polynomial .*)
1.1258 -fun count_neg ([]:mv_poly) = 0
1.1259 - | count_neg ((c,e)::xs) = if c<0 then 1+count_neg xs
1.1260 - else count_neg xs;
1.1261 -
1.1262 -(*. help function for is_polynomial
1.1263 - checks the order of the operators .*)
1.1264 -fun test_polynomial (Const ("uminus",_) $ Free (str,_)) _ = true (*WN.13.3.03*)
1.1265 - | test_polynomial (t as Free(str,_)) v = true
1.1266 - | test_polynomial (t as Const ("op *",_) $ t1 $ t2) v = if v="^" then false
1.1267 - else (test_polynomial t1 "*") andalso (test_polynomial t2 "*")
1.1268 - | test_polynomial (t as Const ("op +",_) $ t1 $ t2) v = if v="*" orelse v="^" then false
1.1269 - else (test_polynomial t1 " ") andalso (test_polynomial t2 " ")
1.1270 - | test_polynomial (t as Const ("Atools.pow",_) $ t1 $ t2) v = (test_polynomial t1 "^") andalso (test_polynomial t2 "^")
1.1271 - | test_polynomial _ v = false;
1.1272 -
1.1273 -(*. tests if a term is a polynomial .*)
1.1274 -fun is_polynomial t = test_polynomial t " ";
1.1275 -
1.1276 -(*. help function for is_expanded
1.1277 - checks the order of the operators .*)
1.1278 -fun test_exp (t as Free(str,_)) v = true
1.1279 - | test_exp (t as Const ("op *",_) $ t1 $ t2) v = if v="^" then false
1.1280 - else (test_exp t1 "*") andalso (test_exp t2 "*")
1.1281 - | test_exp (t as Const ("op +",_) $ t1 $ t2) v = if v="*" orelse v="^" then false
1.1282 - else (test_exp t1 " ") andalso (test_exp t2 " ")
1.1283 - | test_exp (t as Const ("op -",_) $ t1 $ t2) v = if v="*" orelse v="^" then false
1.1284 - else (test_exp t1 " ") andalso (test_exp t2 " ")
1.1285 - | test_exp (t as Const ("Atools.pow",_) $ t1 $ t2) v = (test_exp t1 "^") andalso (test_exp t2 "^")
1.1286 - | test_exp _ v = false;
1.1287 -
1.1288 -
1.1289 -(*. help function for check_coeff:
1.1290 - converts the term to a list of coefficients .*)
1.1291 -fun term2coef' (t as Free(str,_(*typ*))) v :mv_poly option =
1.1292 - let
1.1293 - val x=ref NONE;
1.1294 - val len=ref 0;
1.1295 - val vl=ref [];
1.1296 - val vh=ref [];
1.1297 - val i=ref 0;
1.1298 - in
1.1299 - if is_numeral str then
1.1300 - (
1.1301 - SOME [(((the o int_of_str) str),mv_null2(v))] handle _ => NONE
1.1302 - )
1.1303 - else (* variable *)
1.1304 - (
1.1305 - len:=length(v);
1.1306 - vh:=v;
1.1307 - while ((!len)>(!i)) do
1.1308 - (
1.1309 - if str=hd((!vh)) then
1.1310 - (
1.1311 - vl:=1::(!vl)
1.1312 - )
1.1313 - else
1.1314 - (
1.1315 - vl:=0::(!vl)
1.1316 - );
1.1317 - vh:=tl(!vh);
1.1318 - i:=(!i)+1
1.1319 - );
1.1320 - SOME [(1,rev(!vl))] handle _ => NONE
1.1321 - )
1.1322 - end
1.1323 - | term2coef' (Const ("op *",_) $ t1 $ t2) v :mv_poly option=
1.1324 - let
1.1325 - val t1pp=ref [];
1.1326 - val t2pp=ref [];
1.1327 - val t1c=ref 0;
1.1328 - val t2c=ref 0;
1.1329 - in
1.1330 - (
1.1331 - t1pp:=(#2(hd(the(term2coef' t1 v))));
1.1332 - t2pp:=(#2(hd(the(term2coef' t2 v))));
1.1333 - t1c:=(#1(hd(the(term2coef' t1 v))));
1.1334 - t2c:=(#1(hd(the(term2coef' t2 v))));
1.1335 -
1.1336 - SOME [( (!t1c)*(!t2c) ,( (map op+) ((!t1pp)~~(!t2pp)) ) )] handle _ => NONE
1.1337 -
1.1338 - )
1.1339 - end
1.1340 - | term2coef' (Const ("Atools.pow",_) $ (t1 as Free (str1,_)) $ (t2 as Free (str2,_))) v :mv_poly option=
1.1341 - let
1.1342 - val x=ref NONE;
1.1343 - val len=ref 0;
1.1344 - val vl=ref [];
1.1345 - val vh=ref [];
1.1346 - val vtemp=ref [];
1.1347 - val i=ref 0;
1.1348 - in
1.1349 - (
1.1350 - if (not o is_numeral) str1 andalso is_numeral str2 then
1.1351 - (
1.1352 - len:=length(v);
1.1353 - vh:=v;
1.1354 -
1.1355 - while ((!len)>(!i)) do
1.1356 - (
1.1357 - if str1=hd((!vh)) then
1.1358 - (
1.1359 - vl:=((the o int_of_str) str2)::(!vl)
1.1360 - )
1.1361 - else
1.1362 - (
1.1363 - vl:=0::(!vl)
1.1364 - );
1.1365 - vh:=tl(!vh);
1.1366 - i:=(!i)+1
1.1367 - );
1.1368 - SOME [(1,rev(!vl))] handle _ => NONE
1.1369 - )
1.1370 - else raise error ("RATIONALS_TERM2COEF_EXCEPTION 1: Invalid term")
1.1371 - )
1.1372 - end
1.1373 - | term2coef' (Const ("op +",_) $ t1 $ t2) v :mv_poly option=
1.1374 - (
1.1375 - SOME ((the(term2coef' t1 v)) @ (the(term2coef' t2 v))) handle _ => NONE
1.1376 - )
1.1377 - | term2coef' (Const ("op -",_) $ t1 $ t2) v :mv_poly option=
1.1378 - (
1.1379 - SOME ((the(term2coef' t1 v)) @ mv_skalar_mul((the(term2coef' t2 v)),1)) handle _ => NONE
1.1380 - )
1.1381 - | term2coef' (term) v = raise error ("RATIONALS_TERM2COEF_EXCEPTION 2: Invalid term");
1.1382 -
1.1383 -(*. checks if all coefficients of a polynomial are positiv (except the first) .*)
1.1384 -fun check_coeff t = (* erste Koeffizient kann <0 sein !!! *)
1.1385 - if count_neg(tl(the(term2coef' t (get_vars(t)))))=0 then true
1.1386 - else false;
1.1387 -
1.1388 -(*. checks for expanded term [3] .*)
1.1389 -fun is_expanded t = test_exp t " " andalso check_coeff(t);
1.1390 -
1.1391 -(*WN.7.3.03 Hilfsfunktion f"ur term2poly'*)
1.1392 -fun mk_monom v' p vs =
1.1393 - let fun conv p (v: string) = if v'= v then p else 0
1.1394 - in map (conv p) vs end;
1.1395 -(* mk_monom "y" 5 ["a","b","x","y","z"];
1.1396 -val it = [0,0,0,5,0] : int list*)
1.1397 -
1.1398 -(*. this function converts the term representation into the internal representation mv_poly .*)
1.1399 -fun term2poly' (Const ("uminus",_) $ Free (str,_)) v = (*WN.7.3.03*)
1.1400 - if is_numeral str
1.1401 - then SOME [((the o int_of_str) ("-"^str), mk_monom "#" 0 v)]
1.1402 - else SOME [(~1, mk_monom str 1 v)]
1.1403 -
1.1404 - | term2poly' (Free(str,_)) v :mv_poly option =
1.1405 - let
1.1406 - val x=ref NONE;
1.1407 - val len=ref 0;
1.1408 - val vl=ref [];
1.1409 - val vh=ref [];
1.1410 - val i=ref 0;
1.1411 - in
1.1412 - if is_numeral str then
1.1413 - (
1.1414 - SOME [(((the o int_of_str) str),mv_null2 v)] handle _ => NONE
1.1415 - )
1.1416 - else (* variable *)
1.1417 - (
1.1418 - len:=length v;
1.1419 - vh:= v;
1.1420 - while ((!len)>(!i)) do
1.1421 - (
1.1422 - if str=hd((!vh)) then
1.1423 - (
1.1424 - vl:=1::(!vl)
1.1425 - )
1.1426 - else
1.1427 - (
1.1428 - vl:=0::(!vl)
1.1429 - );
1.1430 - vh:=tl(!vh);
1.1431 - i:=(!i)+1
1.1432 - );
1.1433 - SOME [(1,rev(!vl))] handle _ => NONE
1.1434 - )
1.1435 - end
1.1436 - | term2poly' (Const ("op *",_) $ t1 $ t2) v :mv_poly option=
1.1437 - let
1.1438 - val t1pp=ref [];
1.1439 - val t2pp=ref [];
1.1440 - val t1c=ref 0;
1.1441 - val t2c=ref 0;
1.1442 - in
1.1443 - (
1.1444 - t1pp:=(#2(hd(the(term2poly' t1 v))));
1.1445 - t2pp:=(#2(hd(the(term2poly' t2 v))));
1.1446 - t1c:=(#1(hd(the(term2poly' t1 v))));
1.1447 - t2c:=(#1(hd(the(term2poly' t2 v))));
1.1448 -
1.1449 - SOME [( (!t1c)*(!t2c) ,( (map op+) ((!t1pp)~~(!t2pp)) ) )]
1.1450 - handle _ => NONE
1.1451 -
1.1452 - )
1.1453 - end
1.1454 - | term2poly' (Const ("Atools.pow",_) $ (t1 as Free (str1,_)) $
1.1455 - (t2 as Free (str2,_))) v :mv_poly option=
1.1456 - let
1.1457 - val x=ref NONE;
1.1458 - val len=ref 0;
1.1459 - val vl=ref [];
1.1460 - val vh=ref [];
1.1461 - val vtemp=ref [];
1.1462 - val i=ref 0;
1.1463 - in
1.1464 - (
1.1465 - if (not o is_numeral) str1 andalso is_numeral str2 then
1.1466 - (
1.1467 - len:=length(v);
1.1468 - vh:=v;
1.1469 -
1.1470 - while ((!len)>(!i)) do
1.1471 - (
1.1472 - if str1=hd((!vh)) then
1.1473 - (
1.1474 - vl:=((the o int_of_str) str2)::(!vl)
1.1475 - )
1.1476 - else
1.1477 - (
1.1478 - vl:=0::(!vl)
1.1479 - );
1.1480 - vh:=tl(!vh);
1.1481 - i:=(!i)+1
1.1482 - );
1.1483 - SOME [(1,rev(!vl))] handle _ => NONE
1.1484 - )
1.1485 - else raise error ("RATIONALS_TERM2POLY_EXCEPTION 1: Invalid term")
1.1486 - )
1.1487 - end
1.1488 - | term2poly' (Const ("op +",_) $ t1 $ t2) v :mv_poly option =
1.1489 - (
1.1490 - SOME ((the(term2poly' t1 v)) @ (the(term2poly' t2 v))) handle _ => NONE
1.1491 - )
1.1492 - | term2poly' (Const ("op -",_) $ t1 $ t2) v :mv_poly option =
1.1493 - (
1.1494 - SOME ((the(term2poly' t1 v)) @ mv_skalar_mul((the(term2poly' t2 v)),~1)) handle _ => NONE
1.1495 - )
1.1496 - | term2poly' (term) v = raise error ("RATIONALS_TERM2POLY_EXCEPTION 2: Invalid term");
1.1497 -
1.1498 -(*. translates an Isabelle term into internal representation.
1.1499 - term2poly
1.1500 - fn : term -> (*normalform [2] *)
1.1501 - string list -> (*for ...!!! BITTE DIE ERKLÄRUNG,
1.1502 - DIE DU MIR LETZTES MAL GEGEBEN HAST*)
1.1503 - mv_monom list (*internal representation *)
1.1504 - option (*the translation may fail with NONE*)
1.1505 -.*)
1.1506 -fun term2poly (t:term) v =
1.1507 - if is_polynomial t then term2poly' t v
1.1508 - else raise error ("term2poly: invalid = "^(term2str t));
1.1509 -
1.1510 -(*. same as term2poly with automatic detection of the variables .*)
1.1511 -fun term2polyx t = term2poly t (((map free2str) o vars) t);
1.1512 -
1.1513 -(*. checks if the term is in expanded polynomial form and converts it into the internal representation .*)
1.1514 -fun expanded2poly (t:term) v =
1.1515 - (*if is_expanded t then*) term2poly' t v
1.1516 - (*else raise error ("RATIONALS_EXPANDED2POLY_EXCEPTION: Invalid Polynomial")*);
1.1517 -
1.1518 -(*. same as expanded2poly with automatic detection of the variables .*)
1.1519 -fun expanded2polyx t = expanded2poly t (((map free2str) o vars) t);
1.1520 -
1.1521 -(*. converts a powerproduct into term representation .*)
1.1522 -fun powerproduct2term(xs,v) =
1.1523 - let
1.1524 - val xss=ref xs;
1.1525 - val vv=ref v;
1.1526 - in
1.1527 - (
1.1528 - while hd(!xss)=0 do
1.1529 - (
1.1530 - xss:=tl(!xss);
1.1531 - vv:=tl(!vv)
1.1532 - );
1.1533 -
1.1534 - if list_is_null(tl(!xss)) then
1.1535 - (
1.1536 - if hd(!xss)=1 then Free(hd(!vv), HOLogic.realT)
1.1537 - else
1.1538 - (
1.1539 - Const("Atools.pow",[HOLogic.realT,HOLogic.realT]--->HOLogic.realT) $
1.1540 - Free(hd(!vv), HOLogic.realT) $ Free(str_of_int (hd(!xss)),HOLogic.realT)
1.1541 - )
1.1542 - )
1.1543 - else
1.1544 - (
1.1545 - if hd(!xss)=1 then
1.1546 - (
1.1547 - Const("op *",[HOLogic.realT,HOLogic.realT]--->HOLogic.realT) $
1.1548 - Free(hd(!vv), HOLogic.realT) $
1.1549 - powerproduct2term(tl(!xss),tl(!vv))
1.1550 - )
1.1551 - else
1.1552 - (
1.1553 - Const("op *",[HOLogic.realT,HOLogic.realT]--->HOLogic.realT) $
1.1554 - (
1.1555 - Const("Atools.pow",[HOLogic.realT,HOLogic.realT]--->HOLogic.realT) $
1.1556 - Free(hd(!vv), HOLogic.realT) $ Free(str_of_int (hd(!xss)),HOLogic.realT)
1.1557 - ) $
1.1558 - powerproduct2term(tl(!xss),tl(!vv))
1.1559 - )
1.1560 - )
1.1561 - )
1.1562 - end;
1.1563 -
1.1564 -(*. converts a monom into term representation .*)
1.1565 -(*fun monom2term ((c,e):mv_monom, v:string list) =
1.1566 - if c=0 then Free(str_of_int 0,HOLogic.realT)
1.1567 - else
1.1568 - (
1.1569 - if list_is_null(e) then
1.1570 - (
1.1571 - Free(str_of_int c,HOLogic.realT)
1.1572 - )
1.1573 - else
1.1574 - (
1.1575 - if c=1 then
1.1576 - (
1.1577 - powerproduct2term(e,v)
1.1578 - )
1.1579 - else
1.1580 - (
1.1581 - Const("op *",[HOLogic.realT,HOLogic.realT]--->HOLogic.realT) $
1.1582 - Free(str_of_int c,HOLogic.realT) $
1.1583 - powerproduct2term(e,v)
1.1584 - )
1.1585 - )
1.1586 - );*)
1.1587 -
1.1588 -
1.1589 -(*fun monom2term ((i, is):mv_monom, v) =
1.1590 - if list_is_null is
1.1591 - then
1.1592 - if i >= 0
1.1593 - then Free (str_of_int i, HOLogic.realT)
1.1594 - else Const ("uminus", HOLogic.realT --> HOLogic.realT) $
1.1595 - Free ((str_of_int o abs) i, HOLogic.realT)
1.1596 - else
1.1597 - if i > 0
1.1598 - then Const ("op *", [HOLogic.realT,HOLogic.realT]---> HOLogic.realT) $
1.1599 - (Free (str_of_int i, HOLogic.realT)) $
1.1600 - powerproduct2term(is, v)
1.1601 - else Const ("op *", [HOLogic.realT,HOLogic.realT]---> HOLogic.realT) $
1.1602 - (Const ("uminus", HOLogic.realT --> HOLogic.realT) $
1.1603 - Free ((str_of_int o abs) i, HOLogic.realT)) $
1.1604 - powerproduct2term(is, vs);---------------------------*)
1.1605 -fun monom2term ((i, is) : mv_monom, vs) =
1.1606 - if list_is_null is
1.1607 - then Free (str_of_int i, HOLogic.realT)
1.1608 - else if i = 1
1.1609 - then powerproduct2term (is, vs)
1.1610 - else Const ("op *", [HOLogic.realT, HOLogic.realT] ---> HOLogic.realT) $
1.1611 - (Free (str_of_int i, HOLogic.realT)) $
1.1612 - powerproduct2term (is, vs);
1.1613 -
1.1614 -(*. converts the internal polynomial representation into an Isabelle term.*)
1.1615 -fun poly2term' ([] : mv_poly, vs) = Free(str_of_int 0, HOLogic.realT)
1.1616 - | poly2term' ([(c, e) : mv_monom], vs) = monom2term ((c, e), vs)
1.1617 - | poly2term' ((c, e) :: ces, vs) =
1.1618 - Const("op +", [HOLogic.realT, HOLogic.realT] ---> HOLogic.realT) $
1.1619 - poly2term (ces, vs) $ monom2term ((c, e), vs)
1.1620 -and poly2term (xs, vs) = poly2term' (rev (sort (mv_geq LEX_) (xs)), vs);
1.1621 -
1.1622 -
1.1623 -(*. converts a monom into term representation .*)
1.1624 -(*. ignores the sign of the coefficients => use only for exp-poly functions .*)
1.1625 -fun monom2term2((c,e):mv_monom, v:string list) =
1.1626 - if c=0 then Free(str_of_int 0,HOLogic.realT)
1.1627 - else
1.1628 - (
1.1629 - if list_is_null(e) then
1.1630 - (
1.1631 - Free(str_of_int (abs(c)),HOLogic.realT)
1.1632 - )
1.1633 - else
1.1634 - (
1.1635 - if abs(c)=1 then
1.1636 - (
1.1637 - powerproduct2term(e,v)
1.1638 - )
1.1639 - else
1.1640 - (
1.1641 - Const("op *",[HOLogic.realT,HOLogic.realT]--->HOLogic.realT) $
1.1642 - Free(str_of_int (abs(c)),HOLogic.realT) $
1.1643 - powerproduct2term(e,v)
1.1644 - )
1.1645 - )
1.1646 - );
1.1647 -
1.1648 -(*. converts the expanded polynomial representation into the term representation .*)
1.1649 -fun exp2term' ([]:mv_poly,vars) = Free(str_of_int 0,HOLogic.realT)
1.1650 - | exp2term' ([(c,e)],vars) = monom2term((c,e),vars)
1.1651 - | exp2term' ((c1,e1)::others,vars) =
1.1652 - if c1<0 then
1.1653 - Const("op -",[HOLogic.realT,HOLogic.realT]--->HOLogic.realT) $
1.1654 - exp2term'(others,vars) $
1.1655 - (
1.1656 - monom2term2((c1,e1),vars)
1.1657 - )
1.1658 - else
1.1659 - Const("op +",[HOLogic.realT,HOLogic.realT]--->HOLogic.realT) $
1.1660 - exp2term'(others,vars) $
1.1661 - (
1.1662 - monom2term2((c1,e1),vars)
1.1663 - );
1.1664 -
1.1665 -(*. sorts the powerproduct by lexicographic termorder and converts them into
1.1666 - a term in polynomial representation .*)
1.1667 -fun poly2expanded (xs,vars) = exp2term'(rev(sort (mv_geq LEX_) (xs)),vars);
1.1668 -
1.1669 -(*. converts a polynomial into expanded form .*)
1.1670 -fun polynomial2expanded t =
1.1671 - (let
1.1672 - val vars=(((map free2str) o vars) t);
1.1673 - in
1.1674 - SOME (poly2expanded (the (term2poly t vars), vars))
1.1675 - end) handle _ => NONE;
1.1676 -
1.1677 -(*. converts a polynomial into polynomial form .*)
1.1678 -fun expanded2polynomial t =
1.1679 - (let
1.1680 - val vars=(((map free2str) o vars) t);
1.1681 - in
1.1682 - SOME (poly2term (the (expanded2poly t vars), vars))
1.1683 - end) handle _ => NONE;
1.1684 -
1.1685 -
1.1686 -(*. calculates the greatest common divisor of numerator and denominator and seperates it from each .*)
1.1687 -fun step_cancel (t as Const ("HOL.divide",_) $ p1 $ p2) =
1.1688 - let
1.1689 - val p1' = ref [];
1.1690 - val p2' = ref [];
1.1691 - val p3 = ref []
1.1692 - val vars = rev(get_vars(p1) union get_vars(p2));
1.1693 - in
1.1694 - (
1.1695 - p1':= sort (mv_geq LEX_) (the (term2poly p1 vars ));
1.1696 - p2':= sort (mv_geq LEX_) (the (term2poly p2 vars ));
1.1697 - p3:= sort (mv_geq LEX_) (mv_gcd (!p1') (!p2'));
1.1698 - if (!p3)=[(1,mv_null2(vars))] then
1.1699 - (
1.1700 - Const ("HOL.divide",[HOLogic.realT,HOLogic.realT]--->HOLogic.realT) $ p1 $ p2
1.1701 - )
1.1702 - else
1.1703 - (
1.1704 -
1.1705 - p1':=sort (mv_geq LEX_) (#1(mv_division((!p1'),(!p3),LEX_)));
1.1706 - p2':=sort (mv_geq LEX_) (#1(mv_division((!p2'),(!p3),LEX_)));
1.1707 -
1.1708 - if #1(hd(sort (mv_geq LEX_) (!p2'))) (*mv_lc2(!p2',LEX_)*)>0 then
1.1709 - (
1.1710 - Const ("HOL.divide",[HOLogic.realT,HOLogic.realT]--->HOLogic.realT)
1.1711 - $
1.1712 - (
1.1713 - Const ("op *",[HOLogic.realT,HOLogic.realT]--->HOLogic.realT) $
1.1714 - poly2term(!p1',vars) $
1.1715 - poly2term(!p3,vars)
1.1716 - )
1.1717 - $
1.1718 - (
1.1719 - Const ("op *",[HOLogic.realT,HOLogic.realT]--->HOLogic.realT) $
1.1720 - poly2term(!p2',vars) $
1.1721 - poly2term(!p3,vars)
1.1722 - )
1.1723 - )
1.1724 - else
1.1725 - (
1.1726 - p1':=mv_skalar_mul(!p1',~1);
1.1727 - p2':=mv_skalar_mul(!p2',~1);
1.1728 - p3:=mv_skalar_mul(!p3,~1);
1.1729 - (
1.1730 - Const ("HOL.divide",[HOLogic.realT,HOLogic.realT]--->HOLogic.realT)
1.1731 - $
1.1732 - (
1.1733 - Const ("op *",[HOLogic.realT,HOLogic.realT]--->HOLogic.realT) $
1.1734 - poly2term(!p1',vars) $
1.1735 - poly2term(!p3,vars)
1.1736 - )
1.1737 - $
1.1738 - (
1.1739 - Const ("op *",[HOLogic.realT,HOLogic.realT]--->HOLogic.realT) $
1.1740 - poly2term(!p2',vars) $
1.1741 - poly2term(!p3,vars)
1.1742 - )
1.1743 - )
1.1744 - )
1.1745 - )
1.1746 - )
1.1747 - end
1.1748 -| step_cancel _ = raise error ("RATIONALS_STEP_CANCEL_EXCEPTION: Invalid fraction");
1.1749 -
1.1750 -
1.1751 -(*. same as step_cancel, this time for expanded forms (input+output) .*)
1.1752 -fun step_cancel_expanded (t as Const ("HOL.divide",_) $ p1 $ p2) =
1.1753 - let
1.1754 - val p1' = ref [];
1.1755 - val p2' = ref [];
1.1756 - val p3 = ref []
1.1757 - val vars = rev(get_vars(p1) union get_vars(p2));
1.1758 - in
1.1759 - (
1.1760 - p1':= sort (mv_geq LEX_) (the (expanded2poly p1 vars ));
1.1761 - p2':= sort (mv_geq LEX_) (the (expanded2poly p2 vars ));
1.1762 - p3:= sort (mv_geq LEX_) (mv_gcd (!p1') (!p2'));
1.1763 - if (!p3)=[(1,mv_null2(vars))] then
1.1764 - (
1.1765 - Const ("HOL.divide",[HOLogic.realT,HOLogic.realT]--->HOLogic.realT) $ p1 $ p2
1.1766 - )
1.1767 - else
1.1768 - (
1.1769 -
1.1770 - p1':=sort (mv_geq LEX_) (#1(mv_division((!p1'),(!p3),LEX_)));
1.1771 - p2':=sort (mv_geq LEX_) (#1(mv_division((!p2'),(!p3),LEX_)));
1.1772 -
1.1773 - if #1(hd(sort (mv_geq LEX_) (!p2')))(* mv_lc2(!p2',LEX_)*)>0 then
1.1774 - (
1.1775 - Const ("HOL.divide",[HOLogic.realT,HOLogic.realT]--->HOLogic.realT)
1.1776 - $
1.1777 - (
1.1778 - Const ("op *",[HOLogic.realT,HOLogic.realT]--->HOLogic.realT) $
1.1779 - poly2expanded(!p1',vars) $
1.1780 - poly2expanded(!p3,vars)
1.1781 - )
1.1782 - $
1.1783 - (
1.1784 - Const ("op *",[HOLogic.realT,HOLogic.realT]--->HOLogic.realT) $
1.1785 - poly2expanded(!p2',vars) $
1.1786 - poly2expanded(!p3,vars)
1.1787 - )
1.1788 - )
1.1789 - else
1.1790 - (
1.1791 - p1':=mv_skalar_mul(!p1',~1);
1.1792 - p2':=mv_skalar_mul(!p2',~1);
1.1793 - p3:=mv_skalar_mul(!p3,~1);
1.1794 - (
1.1795 - Const ("HOL.divide",[HOLogic.realT,HOLogic.realT]--->HOLogic.realT)
1.1796 - $
1.1797 - (
1.1798 - Const ("op *",[HOLogic.realT,HOLogic.realT]--->HOLogic.realT) $
1.1799 - poly2expanded(!p1',vars) $
1.1800 - poly2expanded(!p3,vars)
1.1801 - )
1.1802 - $
1.1803 - (
1.1804 - Const ("op *",[HOLogic.realT,HOLogic.realT]--->HOLogic.realT) $
1.1805 - poly2expanded(!p2',vars) $
1.1806 - poly2expanded(!p3,vars)
1.1807 - )
1.1808 - )
1.1809 - )
1.1810 - )
1.1811 - )
1.1812 - end
1.1813 -| step_cancel_expanded _ = raise error ("RATIONALS_STEP_CANCEL_EXCEPTION: Invalid fraction");
1.1814 -
1.1815 -(*. calculates the greatest common divisor of numerator and denominator and divides each through it .*)
1.1816 -fun direct_cancel (t as Const ("HOL.divide",_) $ p1 $ p2) =
1.1817 - let
1.1818 - val p1' = ref [];
1.1819 - val p2' = ref [];
1.1820 - val p3 = ref []
1.1821 - val vars = rev(get_vars(p1) union get_vars(p2));
1.1822 - in
1.1823 - (
1.1824 - p1':=sort (mv_geq LEX_) (mv_shorten((the (term2poly p1 vars )),LEX_));
1.1825 - p2':=sort (mv_geq LEX_) (mv_shorten((the (term2poly p2 vars )),LEX_));
1.1826 - p3 :=sort (mv_geq LEX_) (mv_gcd (!p1') (!p2'));
1.1827 -
1.1828 - if (!p3)=[(1,mv_null2(vars))] then
1.1829 - (
1.1830 - (Const ("HOL.divide",[HOLogic.realT,HOLogic.realT]--->HOLogic.realT) $ p1 $ p2,[])
1.1831 - )
1.1832 - else
1.1833 - (
1.1834 - p1':=sort (mv_geq LEX_) (#1(mv_division((!p1'),(!p3),LEX_)));
1.1835 - p2':=sort (mv_geq LEX_) (#1(mv_division((!p2'),(!p3),LEX_)));
1.1836 - if #1(hd(sort (mv_geq LEX_) (!p2'))) (*mv_lc2(!p2',LEX_)*)>0 then
1.1837 - (
1.1838 - (
1.1839 - Const ("HOL.divide",[HOLogic.realT,HOLogic.realT]--->HOLogic.realT)
1.1840 - $
1.1841 - (
1.1842 - poly2term((!p1'),vars)
1.1843 - )
1.1844 - $
1.1845 - (
1.1846 - poly2term((!p2'),vars)
1.1847 - )
1.1848 - )
1.1849 - ,
1.1850 - if mv_grad(!p3)>0 then
1.1851 - [
1.1852 - (
1.1853 - Const ("Not",[bool]--->bool) $
1.1854 - (
1.1855 - Const("op =",[HOLogic.realT,HOLogic.realT]--->bool) $
1.1856 - poly2term((!p3),vars) $
1.1857 - Free("0",HOLogic.realT)
1.1858 - )
1.1859 - )
1.1860 - ]
1.1861 - else
1.1862 - []
1.1863 - )
1.1864 - else
1.1865 - (
1.1866 - p1':=mv_skalar_mul(!p1',~1);
1.1867 - p2':=mv_skalar_mul(!p2',~1);
1.1868 - if length(!p3)> 2*(count_neg(!p3)) then () else p3 :=mv_skalar_mul(!p3,~1);
1.1869 - (
1.1870 - Const ("HOL.divide",[HOLogic.realT,HOLogic.realT]--->HOLogic.realT)
1.1871 - $
1.1872 - (
1.1873 - poly2term((!p1'),vars)
1.1874 - )
1.1875 - $
1.1876 - (
1.1877 - poly2term((!p2'),vars)
1.1878 - )
1.1879 - ,
1.1880 - if mv_grad(!p3)>0 then
1.1881 - [
1.1882 - (
1.1883 - Const ("Not",[bool]--->bool) $
1.1884 - (
1.1885 - Const("op =",[HOLogic.realT,HOLogic.realT]--->bool) $
1.1886 - poly2term((!p3),vars) $
1.1887 - Free("0",HOLogic.realT)
1.1888 - )
1.1889 - )
1.1890 - ]
1.1891 - else
1.1892 - []
1.1893 - )
1.1894 - )
1.1895 - )
1.1896 - )
1.1897 - end
1.1898 - | direct_cancel _ = raise error ("RATIONALS_DIRECT_CANCEL_EXCEPTION: Invalid fraction");
1.1899 -
1.1900 -(*. same es direct_cancel, this time for expanded forms (input+output).*)
1.1901 -fun direct_cancel_expanded (t as Const ("HOL.divide",_) $ p1 $ p2) =
1.1902 - let
1.1903 - val p1' = ref [];
1.1904 - val p2' = ref [];
1.1905 - val p3 = ref []
1.1906 - val vars = rev(get_vars(p1) union get_vars(p2));
1.1907 - in
1.1908 - (
1.1909 - p1':=sort (mv_geq LEX_) (mv_shorten((the (expanded2poly p1 vars )),LEX_));
1.1910 - p2':=sort (mv_geq LEX_) (mv_shorten((the (expanded2poly p2 vars )),LEX_));
1.1911 - p3 :=sort (mv_geq LEX_) (mv_gcd (!p1') (!p2'));
1.1912 -
1.1913 - if (!p3)=[(1,mv_null2(vars))] then
1.1914 - (
1.1915 - (Const ("HOL.divide",[HOLogic.realT,HOLogic.realT]--->HOLogic.realT) $ p1 $ p2,[])
1.1916 - )
1.1917 - else
1.1918 - (
1.1919 - p1':=sort (mv_geq LEX_) (#1(mv_division((!p1'),(!p3),LEX_)));
1.1920 - p2':=sort (mv_geq LEX_) (#1(mv_division((!p2'),(!p3),LEX_)));
1.1921 - if #1(hd(sort (mv_geq LEX_) (!p2'))) (*mv_lc2(!p2',LEX_)*)>0 then
1.1922 - (
1.1923 - (
1.1924 - Const ("HOL.divide",[HOLogic.realT,HOLogic.realT]--->HOLogic.realT)
1.1925 - $
1.1926 - (
1.1927 - poly2expanded((!p1'),vars)
1.1928 - )
1.1929 - $
1.1930 - (
1.1931 - poly2expanded((!p2'),vars)
1.1932 - )
1.1933 - )
1.1934 - ,
1.1935 - if mv_grad(!p3)>0 then
1.1936 - [
1.1937 - (
1.1938 - Const ("Not",[bool]--->bool) $
1.1939 - (
1.1940 - Const("op =",[HOLogic.realT,HOLogic.realT]--->bool) $
1.1941 - poly2expanded((!p3),vars) $
1.1942 - Free("0",HOLogic.realT)
1.1943 - )
1.1944 - )
1.1945 - ]
1.1946 - else
1.1947 - []
1.1948 - )
1.1949 - else
1.1950 - (
1.1951 - p1':=mv_skalar_mul(!p1',~1);
1.1952 - p2':=mv_skalar_mul(!p2',~1);
1.1953 - if length(!p3)> 2*(count_neg(!p3)) then () else p3 :=mv_skalar_mul(!p3,~1);
1.1954 - (
1.1955 - Const ("HOL.divide",[HOLogic.realT,HOLogic.realT]--->HOLogic.realT)
1.1956 - $
1.1957 - (
1.1958 - poly2expanded((!p1'),vars)
1.1959 - )
1.1960 - $
1.1961 - (
1.1962 - poly2expanded((!p2'),vars)
1.1963 - )
1.1964 - ,
1.1965 - if mv_grad(!p3)>0 then
1.1966 - [
1.1967 - (
1.1968 - Const ("Not",[bool]--->bool) $
1.1969 - (
1.1970 - Const("op =",[HOLogic.realT,HOLogic.realT]--->bool) $
1.1971 - poly2expanded((!p3),vars) $
1.1972 - Free("0",HOLogic.realT)
1.1973 - )
1.1974 - )
1.1975 - ]
1.1976 - else
1.1977 - []
1.1978 - )
1.1979 - )
1.1980 - )
1.1981 - )
1.1982 - end
1.1983 - | direct_cancel_expanded _ = raise error ("RATIONALS_DIRECT_CANCEL_EXCEPTION: Invalid fraction");
1.1984 -
1.1985 -
1.1986 -(*. adds two fractions .*)
1.1987 -fun add_fract ((Const("HOL.divide",_) $ t11 $ t12),(Const("HOL.divide",_) $ t21 $ t22)) =
1.1988 - let
1.1989 - val vars=get_vars(t11) union get_vars(t12) union get_vars(t21) union get_vars(t22);
1.1990 - val t11'=ref (the(term2poly t11 vars));
1.1991 -val _= writeln"### add_fract: done t11"
1.1992 - val t12'=ref (the(term2poly t12 vars));
1.1993 -val _= writeln"### add_fract: done t12"
1.1994 - val t21'=ref (the(term2poly t21 vars));
1.1995 -val _= writeln"### add_fract: done t21"
1.1996 - val t22'=ref (the(term2poly t22 vars));
1.1997 -val _= writeln"### add_fract: done t22"
1.1998 - val den=ref [];
1.1999 - val nom=ref [];
1.2000 - val m1=ref [];
1.2001 - val m2=ref [];
1.2002 - in
1.2003 -
1.2004 - (
1.2005 - den :=sort (mv_geq LEX_) (mv_lcm (!t12') (!t22'));
1.2006 -writeln"### add_fract: done sort mv_lcm";
1.2007 - m1 :=sort (mv_geq LEX_) (#1(mv_division(!den,!t12',LEX_)));
1.2008 -writeln"### add_fract: done sort mv_division t12";
1.2009 - m2 :=sort (mv_geq LEX_) (#1(mv_division(!den,!t22',LEX_)));
1.2010 -writeln"### add_fract: done sort mv_division t22";
1.2011 - nom :=sort (mv_geq LEX_)
1.2012 - (mv_shorten(mv_add(mv_mul(!t11',!m1,LEX_),
1.2013 - mv_mul(!t21',!m2,LEX_),
1.2014 - LEX_),
1.2015 - LEX_));
1.2016 -writeln"### add_fract: done sort mv_add";
1.2017 - (
1.2018 - Const ("HOL.divide",[HOLogic.realT,HOLogic.realT]--->HOLogic.realT)
1.2019 - $
1.2020 - (
1.2021 - poly2term((!nom),vars)
1.2022 - )
1.2023 - $
1.2024 - (
1.2025 - poly2term((!den),vars)
1.2026 - )
1.2027 - )
1.2028 - )
1.2029 - end
1.2030 - | add_fract (_,_) = raise error ("RATIONALS_ADD_FRACTION_EXCEPTION: Invalid add_fraction call");
1.2031 -
1.2032 -(*. adds two expanded fractions .*)
1.2033 -fun add_fract_exp ((Const("HOL.divide",_) $ t11 $ t12),(Const("HOL.divide",_) $ t21 $ t22)) =
1.2034 - let
1.2035 - val vars=get_vars(t11) union get_vars(t12) union get_vars(t21) union get_vars(t22);
1.2036 - val t11'=ref (the(expanded2poly t11 vars));
1.2037 - val t12'=ref (the(expanded2poly t12 vars));
1.2038 - val t21'=ref (the(expanded2poly t21 vars));
1.2039 - val t22'=ref (the(expanded2poly t22 vars));
1.2040 - val den=ref [];
1.2041 - val nom=ref [];
1.2042 - val m1=ref [];
1.2043 - val m2=ref [];
1.2044 - in
1.2045 -
1.2046 - (
1.2047 - den :=sort (mv_geq LEX_) (mv_lcm (!t12') (!t22'));
1.2048 - m1 :=sort (mv_geq LEX_) (#1(mv_division(!den,!t12',LEX_)));
1.2049 - m2 :=sort (mv_geq LEX_) (#1(mv_division(!den,!t22',LEX_)));
1.2050 - nom :=sort (mv_geq LEX_) (mv_shorten(mv_add(mv_mul(!t11',!m1,LEX_),mv_mul(!t21',!m2,LEX_),LEX_),LEX_));
1.2051 - (
1.2052 - Const ("HOL.divide",[HOLogic.realT,HOLogic.realT]--->HOLogic.realT)
1.2053 - $
1.2054 - (
1.2055 - poly2expanded((!nom),vars)
1.2056 - )
1.2057 - $
1.2058 - (
1.2059 - poly2expanded((!den),vars)
1.2060 - )
1.2061 - )
1.2062 - )
1.2063 - end
1.2064 - | add_fract_exp (_,_) = raise error ("RATIONALS_ADD_FRACTION_EXP_EXCEPTION: Invalid add_fraction call");
1.2065 -
1.2066 -(*. adds a list of terms .*)
1.2067 -fun add_list_of_fractions []= (Free("0",HOLogic.realT),[])
1.2068 - | add_list_of_fractions [x]= direct_cancel x
1.2069 - | add_list_of_fractions (x::y::xs) =
1.2070 - let
1.2071 - val (t1a,rest1)=direct_cancel(x);
1.2072 -val _= writeln"### add_list_of_fractions xs: has done direct_cancel(x)";
1.2073 - val (t2a,rest2)=direct_cancel(y);
1.2074 -val _= writeln"### add_list_of_fractions xs: has done direct_cancel(y)";
1.2075 - val (t3a,rest3)=(add_list_of_fractions (add_fract(t1a,t2a)::xs));
1.2076 -val _= writeln"### add_list_of_fractions xs: has done add_list_of_fraction xs";
1.2077 - val (t4a,rest4)=direct_cancel(t3a);
1.2078 -val _= writeln"### add_list_of_fractions xs: has done direct_cancel(t3a)";
1.2079 - val rest=rest1 union rest2 union rest3 union rest4;
1.2080 - in
1.2081 - (writeln"### add_list_of_fractions in";
1.2082 - (
1.2083 - (t4a,rest)
1.2084 - )
1.2085 - )
1.2086 - end;
1.2087 -
1.2088 -(*. adds a list of expanded terms .*)
1.2089 -fun add_list_of_fractions_exp []= (Free("0",HOLogic.realT),[])
1.2090 - | add_list_of_fractions_exp [x]= direct_cancel_expanded x
1.2091 - | add_list_of_fractions_exp (x::y::xs) =
1.2092 - let
1.2093 - val (t1a,rest1)=direct_cancel_expanded(x);
1.2094 - val (t2a,rest2)=direct_cancel_expanded(y);
1.2095 - val (t3a,rest3)=(add_list_of_fractions_exp (add_fract_exp(t1a,t2a)::xs));
1.2096 - val (t4a,rest4)=direct_cancel_expanded(t3a);
1.2097 - val rest=rest1 union rest2 union rest3 union rest4;
1.2098 - in
1.2099 - (
1.2100 - (t4a,rest)
1.2101 - )
1.2102 - end;
1.2103 -
1.2104 -(*. calculates the lcm of a list of mv_poly .*)
1.2105 -fun calc_lcm ([x],var)= (x,var)
1.2106 - | calc_lcm ((x::xs),var) = (mv_lcm x (#1(calc_lcm (xs,var))),var);
1.2107 -
1.2108 -(*. converts a list of terms to a list of mv_poly .*)
1.2109 -fun t2d([],_)=[]
1.2110 - | t2d((t as (Const("HOL.divide",_) $ p1 $ p2))::xs,vars)= (the(term2poly p2 vars)) :: t2d(xs,vars);
1.2111 -
1.2112 -(*. same as t2d, this time for expanded forms .*)
1.2113 -fun t2d_exp([],_)=[]
1.2114 - | t2d_exp((t as (Const("HOL.divide",_) $ p1 $ p2))::xs,vars)= (the(expanded2poly p2 vars)) :: t2d_exp(xs,vars);
1.2115 -
1.2116 -(*. converts a list of fract terms to a list of their denominators .*)
1.2117 -fun termlist2denominators [] = ([],[])
1.2118 - | termlist2denominators xs =
1.2119 - let
1.2120 - val xxs=ref xs;
1.2121 - val var=ref [];
1.2122 - in
1.2123 - var:=[];
1.2124 - while length(!xxs)>0 do
1.2125 - (
1.2126 - let
1.2127 - val (t as Const ("HOL.divide",_) $ p1x $ p2x)=hd(!xxs);
1.2128 - in
1.2129 - (
1.2130 - xxs:=tl(!xxs);
1.2131 - var:=((get_vars(p2x)) union (get_vars(p1x)) union (!var))
1.2132 - )
1.2133 - end
1.2134 - );
1.2135 - (t2d(xs,!var),!var)
1.2136 - end;
1.2137 -
1.2138 -(*. calculates the lcm of a list of mv_poly .*)
1.2139 -fun calc_lcm ([x],var)= (x,var)
1.2140 - | calc_lcm ((x::xs),var) = (mv_lcm x (#1(calc_lcm (xs,var))),var);
1.2141 -
1.2142 -(*. converts a list of terms to a list of mv_poly .*)
1.2143 -fun t2d([],_)=[]
1.2144 - | t2d((t as (Const("HOL.divide",_) $ p1 $ p2))::xs,vars)= (the(term2poly p2 vars)) :: t2d(xs,vars);
1.2145 -
1.2146 -(*. same as t2d, this time for expanded forms .*)
1.2147 -fun t2d_exp([],_)=[]
1.2148 - | t2d_exp((t as (Const("HOL.divide",_) $ p1 $ p2))::xs,vars)= (the(expanded2poly p2 vars)) :: t2d_exp(xs,vars);
1.2149 -
1.2150 -(*. converts a list of fract terms to a list of their denominators .*)
1.2151 -fun termlist2denominators [] = ([],[])
1.2152 - | termlist2denominators xs =
1.2153 - let
1.2154 - val xxs=ref xs;
1.2155 - val var=ref [];
1.2156 - in
1.2157 - var:=[];
1.2158 - while length(!xxs)>0 do
1.2159 - (
1.2160 - let
1.2161 - val (t as Const ("HOL.divide",_) $ p1x $ p2x)=hd(!xxs);
1.2162 - in
1.2163 - (
1.2164 - xxs:=tl(!xxs);
1.2165 - var:=((get_vars(p2x)) union (get_vars(p1x)) union (!var))
1.2166 - )
1.2167 - end
1.2168 - );
1.2169 - (t2d(xs,!var),!var)
1.2170 - end;
1.2171 -
1.2172 -(*. same as termlist2denminators, this time for expanded forms .*)
1.2173 -fun termlist2denominators_exp [] = ([],[])
1.2174 - | termlist2denominators_exp xs =
1.2175 - let
1.2176 - val xxs=ref xs;
1.2177 - val var=ref [];
1.2178 - in
1.2179 - var:=[];
1.2180 - while length(!xxs)>0 do
1.2181 - (
1.2182 - let
1.2183 - val (t as Const ("HOL.divide",_) $ p1x $ p2x)=hd(!xxs);
1.2184 - in
1.2185 - (
1.2186 - xxs:=tl(!xxs);
1.2187 - var:=((get_vars(p2x)) union (get_vars(p1x)) union (!var))
1.2188 - )
1.2189 - end
1.2190 - );
1.2191 - (t2d_exp(xs,!var),!var)
1.2192 - end;
1.2193 -
1.2194 -(*. reduces all fractions to the least common denominator .*)
1.2195 -fun com_den(x::xs,denom,den,var)=
1.2196 - let
1.2197 - val (t as Const ("HOL.divide",_) $ p1' $ p2')=x;
1.2198 - val p2= sort (mv_geq LEX_) (the(term2poly p2' var));
1.2199 - val p3= #1(mv_division(denom,p2,LEX_));
1.2200 - val p1var=get_vars(p1');
1.2201 - in
1.2202 - if length(xs)>0 then
1.2203 - if p3=[(1,mv_null2(var))] then
1.2204 - (
1.2205 - Const ("op +",[HOLogic.realT,HOLogic.realT]--->HOLogic.realT)
1.2206 - $
1.2207 - (
1.2208 - Const ("HOL.divide",[HOLogic.realT,HOLogic.realT]--->HOLogic.realT)
1.2209 - $
1.2210 - poly2term(the (term2poly p1' p1var),p1var)
1.2211 - $
1.2212 - den
1.2213 - )
1.2214 - $
1.2215 - #1(com_den(xs,denom,den,var))
1.2216 - ,
1.2217 - []
1.2218 - )
1.2219 - else
1.2220 - (
1.2221 - Const ("op +",[HOLogic.realT,HOLogic.realT]--->HOLogic.realT)
1.2222 - $
1.2223 - (
1.2224 - Const ("HOL.divide",[HOLogic.realT,HOLogic.realT]--->HOLogic.realT)
1.2225 - $
1.2226 - (
1.2227 - Const ("op *",[HOLogic.realT,HOLogic.realT]--->HOLogic.realT) $
1.2228 - poly2term(the (term2poly p1' p1var),p1var) $
1.2229 - poly2term(p3,var)
1.2230 - )
1.2231 - $
1.2232 - (
1.2233 - den
1.2234 - )
1.2235 - )
1.2236 - $
1.2237 - #1(com_den(xs,denom,den,var))
1.2238 - ,
1.2239 - []
1.2240 - )
1.2241 - else
1.2242 - if p3=[(1,mv_null2(var))] then
1.2243 - (
1.2244 - (
1.2245 - Const ("HOL.divide",[HOLogic.realT,HOLogic.realT]--->HOLogic.realT)
1.2246 - $
1.2247 - poly2term(the (term2poly p1' p1var),p1var)
1.2248 - $
1.2249 - den
1.2250 - )
1.2251 - ,
1.2252 - []
1.2253 - )
1.2254 - else
1.2255 - (
1.2256 - Const ("HOL.divide",[HOLogic.realT,HOLogic.realT]--->HOLogic.realT)
1.2257 - $
1.2258 - (
1.2259 - Const ("op *",[HOLogic.realT,HOLogic.realT]--->HOLogic.realT) $
1.2260 - poly2term(the (term2poly p1' p1var),p1var) $
1.2261 - poly2term(p3,var)
1.2262 - )
1.2263 - $
1.2264 - den
1.2265 - ,
1.2266 - []
1.2267 - )
1.2268 - end;
1.2269 -
1.2270 -(*. same as com_den, this time for expanded forms .*)
1.2271 -fun com_den_exp(x::xs,denom,den,var)=
1.2272 - let
1.2273 - val (t as Const ("HOL.divide",_) $ p1' $ p2')=x;
1.2274 - val p2= sort (mv_geq LEX_) (the(expanded2poly p2' var));
1.2275 - val p3= #1(mv_division(denom,p2,LEX_));
1.2276 - val p1var=get_vars(p1');
1.2277 - in
1.2278 - if length(xs)>0 then
1.2279 - if p3=[(1,mv_null2(var))] then
1.2280 - (
1.2281 - Const ("op +",[HOLogic.realT,HOLogic.realT]--->HOLogic.realT)
1.2282 - $
1.2283 - (
1.2284 - Const ("HOL.divide",[HOLogic.realT,HOLogic.realT]--->HOLogic.realT)
1.2285 - $
1.2286 - poly2expanded(the(expanded2poly p1' p1var),p1var)
1.2287 - $
1.2288 - den
1.2289 - )
1.2290 - $
1.2291 - #1(com_den_exp(xs,denom,den,var))
1.2292 - ,
1.2293 - []
1.2294 - )
1.2295 - else
1.2296 - (
1.2297 - Const ("op +",[HOLogic.realT,HOLogic.realT]--->HOLogic.realT)
1.2298 - $
1.2299 - (
1.2300 - Const ("HOL.divide",[HOLogic.realT,HOLogic.realT]--->HOLogic.realT)
1.2301 - $
1.2302 - (
1.2303 - Const ("op *",[HOLogic.realT,HOLogic.realT]--->HOLogic.realT) $
1.2304 - poly2expanded(the(expanded2poly p1' p1var),p1var) $
1.2305 - poly2expanded(p3,var)
1.2306 - )
1.2307 - $
1.2308 - (
1.2309 - den
1.2310 - )
1.2311 - )
1.2312 - $
1.2313 - #1(com_den_exp(xs,denom,den,var))
1.2314 - ,
1.2315 - []
1.2316 - )
1.2317 - else
1.2318 - if p3=[(1,mv_null2(var))] then
1.2319 - (
1.2320 - (
1.2321 - Const ("HOL.divide",[HOLogic.realT,HOLogic.realT]--->HOLogic.realT)
1.2322 - $
1.2323 - poly2expanded(the(expanded2poly p1' p1var),p1var)
1.2324 - $
1.2325 - den
1.2326 - )
1.2327 - ,
1.2328 - []
1.2329 - )
1.2330 - else
1.2331 - (
1.2332 - Const ("HOL.divide",[HOLogic.realT,HOLogic.realT]--->HOLogic.realT)
1.2333 - $
1.2334 - (
1.2335 - Const ("op *",[HOLogic.realT,HOLogic.realT]--->HOLogic.realT) $
1.2336 - poly2expanded(the(expanded2poly p1' p1var),p1var) $
1.2337 - poly2expanded(p3,var)
1.2338 - )
1.2339 - $
1.2340 - den
1.2341 - ,
1.2342 - []
1.2343 - )
1.2344 - end;
1.2345 -
1.2346 -(* wird aktuell nicht mehr gebraucht, bei rückänderung schon
1.2347 --------------------------------------------------------------
1.2348 -(* WN0210???SK brauch ma des überhaupt *)
1.2349 -fun com_den2(x::xs,denom,den,var)=
1.2350 - let
1.2351 - val (t as Const ("HOL.divide",_) $ p1' $ p2')=x;
1.2352 - val p2= sort (mv_geq LEX_) (the(term2poly p2' var));
1.2353 - val p3= #1(mv_division(denom,p2,LEX_));
1.2354 - val p1var=get_vars(p1');
1.2355 - in
1.2356 - if length(xs)>0 then
1.2357 - if p3=[(1,mv_null2(var))] then
1.2358 - (
1.2359 - Const ("op +",[HOLogic.realT,HOLogic.realT]--->HOLogic.realT) $
1.2360 - poly2term(the(term2poly p1' p1var),p1var) $
1.2361 - com_den2(xs,denom,den,var)
1.2362 - )
1.2363 - else
1.2364 - (
1.2365 - Const ("op +",[HOLogic.realT,HOLogic.realT]--->HOLogic.realT) $
1.2366 - (
1.2367 - let
1.2368 - val p3'=poly2term(p3,var);
1.2369 - val vars= (((map free2str) o vars) p1') union (((map free2str) o vars) p3');
1.2370 - in
1.2371 - poly2term(sort (mv_geq LEX_) (mv_mul(the(term2poly p1' vars) ,the(term2poly p3' vars),LEX_)),vars)
1.2372 - end
1.2373 - ) $
1.2374 - com_den2(xs,denom,den,var)
1.2375 - )
1.2376 - else
1.2377 - if p3=[(1,mv_null2(var))] then
1.2378 - (
1.2379 - poly2term(the(term2poly p1' p1var),p1var)
1.2380 - )
1.2381 - else
1.2382 - (
1.2383 - let
1.2384 - val p3'=poly2term(p3,var);
1.2385 - val vars= (((map free2str) o vars) p1') union (((map free2str) o vars) p3');
1.2386 - in
1.2387 - poly2term(sort (mv_geq LEX_) (mv_mul(the(term2poly p1' vars) ,the(term2poly p3' vars),LEX_)),vars)
1.2388 - end
1.2389 - )
1.2390 - end;
1.2391 -
1.2392 -(* WN0210???SK brauch ma des überhaupt *)
1.2393 -fun com_den_exp2(x::xs,denom,den,var)=
1.2394 - let
1.2395 - val (t as Const ("HOL.divide",_) $ p1' $ p2')=x;
1.2396 - val p2= sort (mv_geq LEX_) (the(expanded2poly p2' var));
1.2397 - val p3= #1(mv_division(denom,p2,LEX_));
1.2398 - val p1var=get_vars p1';
1.2399 - in
1.2400 - if length(xs)>0 then
1.2401 - if p3=[(1,mv_null2(var))] then
1.2402 - (
1.2403 - Const ("op +",[HOLogic.realT,HOLogic.realT]--->HOLogic.realT) $
1.2404 - poly2expanded(the (expanded2poly p1' p1var),p1var) $
1.2405 - com_den_exp2(xs,denom,den,var)
1.2406 - )
1.2407 - else
1.2408 - (
1.2409 - Const ("op +",[HOLogic.realT,HOLogic.realT]--->HOLogic.realT) $
1.2410 - (
1.2411 - let
1.2412 - val p3'=poly2expanded(p3,var);
1.2413 - val vars= (((map free2str) o vars) p1') union (((map free2str) o vars) p3');
1.2414 - in
1.2415 - poly2expanded(sort (mv_geq LEX_) (mv_mul(the(expanded2poly p1' vars) ,the(expanded2poly p3' vars),LEX_)),vars)
1.2416 - end
1.2417 - ) $
1.2418 - com_den_exp2(xs,denom,den,var)
1.2419 - )
1.2420 - else
1.2421 - if p3=[(1,mv_null2(var))] then
1.2422 - (
1.2423 - poly2expanded(the (expanded2poly p1' p1var),p1var)
1.2424 - )
1.2425 - else
1.2426 - (
1.2427 - let
1.2428 - val p3'=poly2expanded(p3,var);
1.2429 - val vars= (((map free2str) o vars) p1') union (((map free2str) o vars) p3');
1.2430 - in
1.2431 - poly2expanded(sort (mv_geq LEX_) (mv_mul(the(expanded2poly p1' vars) ,the(expanded2poly p3' vars),LEX_)),vars)
1.2432 - end
1.2433 - )
1.2434 - end;
1.2435 ----------------------------------------------------------*)
1.2436 -
1.2437 -
1.2438 -(*. searches for an element y of a list ys, which has an gcd not 1 with x .*)
1.2439 -fun exists_gcd (x,[]) = false
1.2440 - | exists_gcd (x,y::ys) = if mv_gcd x y = [(1,mv_null2(#2(hd(x))))] then exists_gcd (x,ys)
1.2441 - else true;
1.2442 -
1.2443 -(*. divides each element of the list xs with y .*)
1.2444 -fun list_div ([],y) = []
1.2445 - | list_div (x::xs,y) =
1.2446 - let
1.2447 - val (d,r)=mv_division(x,y,LEX_);
1.2448 - in
1.2449 - if r=[] then
1.2450 - d::list_div(xs,y)
1.2451 - else x::list_div(xs,y)
1.2452 - end;
1.2453 -
1.2454 -(*. checks if x is in the list ys .*)
1.2455 -fun in_list (x,[]) = false
1.2456 - | in_list (x,y::ys) = if x=y then true
1.2457 - else in_list(x,ys);
1.2458 -
1.2459 -(*. deletes all equal elements of the list xs .*)
1.2460 -fun kill_equal [] = []
1.2461 - | kill_equal (x::xs) = if in_list(x,xs) orelse x=[(1,mv_null2(#2(hd(x))))] then kill_equal(xs)
1.2462 - else x::kill_equal(xs);
1.2463 -
1.2464 -(*. searches for new factors .*)
1.2465 -fun new_factors [] = []
1.2466 - | new_factors (list:mv_poly list):mv_poly list =
1.2467 - let
1.2468 - val l = kill_equal list;
1.2469 - val len = length(l);
1.2470 - in
1.2471 - if len>=2 then
1.2472 - (
1.2473 - let
1.2474 - val x::y::xs=l;
1.2475 - val gcd=mv_gcd x y;
1.2476 - in
1.2477 - if gcd=[(1,mv_null2(#2(hd(x))))] then
1.2478 - (
1.2479 - if exists_gcd(x,xs) then new_factors (y::xs @ [x])
1.2480 - else x::new_factors(y::xs)
1.2481 - )
1.2482 - else gcd::new_factors(kill_equal(list_div(x::y::xs,gcd)))
1.2483 - end
1.2484 - )
1.2485 - else
1.2486 - if len=1 then [hd(l)]
1.2487 - else []
1.2488 - end;
1.2489 -
1.2490 -(*. gets the factors of a list .*)
1.2491 -fun get_factors x = new_factors x;
1.2492 -
1.2493 -(*. multiplies the elements of the list .*)
1.2494 -fun multi_list [] = []
1.2495 - | multi_list (x::xs) = if xs=[] then x
1.2496 - else mv_mul(x,multi_list xs,LEX_);
1.2497 -
1.2498 -(*. makes a term out of the elements of the list (polynomial representation) .*)
1.2499 -fun make_term ([],vars) = Free(str_of_int 0,HOLogic.realT)
1.2500 - | make_term ((x::xs),vars) = if length(xs)=0 then poly2term(sort (mv_geq LEX_) (x),vars)
1.2501 - else
1.2502 - (
1.2503 - Const ("op *",[HOLogic.realT,HOLogic.realT]--->HOLogic.realT) $
1.2504 - poly2term(sort (mv_geq LEX_) (x),vars) $
1.2505 - make_term(xs,vars)
1.2506 - );
1.2507 -
1.2508 -(*. factorizes the denominator (polynomial representation) .*)
1.2509 -fun factorize_den (l,den,vars) =
1.2510 - let
1.2511 - val factor_list=kill_equal( (get_factors l));
1.2512 - val mlist=multi_list(factor_list);
1.2513 - val (last,rest)=mv_division(den,multi_list(factor_list),LEX_);
1.2514 - in
1.2515 - if rest=[] then
1.2516 - (
1.2517 - if last=[(1,mv_null2(vars))] then make_term(factor_list,vars)
1.2518 - else make_term(last::factor_list,vars)
1.2519 - )
1.2520 - else raise error ("RATIONALS_FACTORIZE_DEN_EXCEPTION: Invalid factor by division")
1.2521 - end;
1.2522 -
1.2523 -(*. makes a term out of the elements of the list (expanded polynomial representation) .*)
1.2524 -fun make_exp ([],vars) = Free(str_of_int 0,HOLogic.realT)
1.2525 - | make_exp ((x::xs),vars) = if length(xs)=0 then poly2expanded(sort (mv_geq LEX_) (x),vars)
1.2526 - else
1.2527 - (
1.2528 - Const ("op *",[HOLogic.realT,HOLogic.realT]--->HOLogic.realT) $
1.2529 - poly2expanded(sort (mv_geq LEX_) (x),vars) $
1.2530 - make_exp(xs,vars)
1.2531 - );
1.2532 -
1.2533 -(*. factorizes the denominator (expanded polynomial representation) .*)
1.2534 -fun factorize_den_exp (l,den,vars) =
1.2535 - let
1.2536 - val factor_list=kill_equal( (get_factors l));
1.2537 - val mlist=multi_list(factor_list);
1.2538 - val (last,rest)=mv_division(den,multi_list(factor_list),LEX_);
1.2539 - in
1.2540 - if rest=[] then
1.2541 - (
1.2542 - if last=[(1,mv_null2(vars))] then make_exp(factor_list,vars)
1.2543 - else make_exp(last::factor_list,vars)
1.2544 - )
1.2545 - else raise error ("RATIONALS_FACTORIZE_DEN_EXP_EXCEPTION: Invalid factor by division")
1.2546 - end;
1.2547 -
1.2548 -(*. calculates the common denominator of all elements of the list and multiplies .*)
1.2549 -(*. the nominators and denominators with the correct factor .*)
1.2550 -(*. (polynomial representation) .*)
1.2551 -fun step_add_list_of_fractions []=(Free("0",HOLogic.realT),[]:term list)
1.2552 - | step_add_list_of_fractions [x]= raise error ("RATIONALS_STEP_ADD_LIST_OF_FRACTIONS_EXCEPTION: Nothing to add")
1.2553 - | step_add_list_of_fractions (xs) =
1.2554 - let
1.2555 - val den_list=termlist2denominators (xs); (* list of denominators *)
1.2556 - val (denom,var)=calc_lcm(den_list); (* common denominator *)
1.2557 - val den=factorize_den(#1(den_list),denom,var); (* faktorisierter Nenner !!! *)
1.2558 - in
1.2559 - com_den(xs,denom,den,var)
1.2560 - end;
1.2561 -
1.2562 -(*. calculates the common denominator of all elements of the list and multiplies .*)
1.2563 -(*. the nominators and denominators with the correct factor .*)
1.2564 -(*. (expanded polynomial representation) .*)
1.2565 -fun step_add_list_of_fractions_exp [] = (Free("0",HOLogic.realT),[]:term list)
1.2566 - | step_add_list_of_fractions_exp [x] = raise error ("RATIONALS_STEP_ADD_LIST_OF_FRACTIONS_EXP_EXCEPTION: Nothing to add")
1.2567 - | step_add_list_of_fractions_exp (xs)=
1.2568 - let
1.2569 - val den_list=termlist2denominators_exp (xs); (* list of denominators *)
1.2570 - val (denom,var)=calc_lcm(den_list); (* common denominator *)
1.2571 - val den=factorize_den_exp(#1(den_list),denom,var); (* faktorisierter Nenner !!! *)
1.2572 - in
1.2573 - com_den_exp(xs,denom,den,var)
1.2574 - end;
1.2575 -
1.2576 -(* wird aktuell nicht mehr gebraucht, bei rückänderung schon
1.2577 --------------------------------------------------------------
1.2578 -(* WN0210???SK brauch ma des überhaupt *)
1.2579 -fun step_add_list_of_fractions2 []=(Free("0",HOLogic.realT),[]:term list)
1.2580 - | step_add_list_of_fractions2 [x]=(x,[])
1.2581 - | step_add_list_of_fractions2 (xs) =
1.2582 - let
1.2583 - val den_list=termlist2denominators (xs); (* list of denominators *)
1.2584 - val (denom,var)=calc_lcm(den_list); (* common denominator *)
1.2585 - val den=factorize_den(#1(den_list),denom,var); (* faktorisierter Nenner !!! *)
1.2586 - in
1.2587 - (
1.2588 - Const ("HOL.divide",[HOLogic.realT,HOLogic.realT]--->HOLogic.realT) $
1.2589 - com_den2(xs,denom, poly2term(denom,var)(*den*),var) $
1.2590 - poly2term(denom,var)
1.2591 - ,
1.2592 - []
1.2593 - )
1.2594 - end;
1.2595 -
1.2596 -(* WN0210???SK brauch ma des überhaupt *)
1.2597 -fun step_add_list_of_fractions2_exp []=(Free("0",HOLogic.realT),[]:term list)
1.2598 - | step_add_list_of_fractions2_exp [x]=(x,[])
1.2599 - | step_add_list_of_fractions2_exp (xs) =
1.2600 - let
1.2601 - val den_list=termlist2denominators_exp (xs); (* list of denominators *)
1.2602 - val (denom,var)=calc_lcm(den_list); (* common denominator *)
1.2603 - val den=factorize_den_exp(#1(den_list),denom,var); (* faktorisierter Nenner !!! *)
1.2604 - in
1.2605 - (
1.2606 - Const ("HOL.divide",[HOLogic.realT,HOLogic.realT]--->HOLogic.realT) $
1.2607 - com_den_exp2(xs,denom, poly2term(denom,var)(*den*),var) $
1.2608 - poly2expanded(denom,var)
1.2609 - ,
1.2610 - []
1.2611 - )
1.2612 - end;
1.2613 ----------------------------------------------- *)
1.2614 -
1.2615 -
1.2616 -(*. converts a term, which contains severel terms seperated by +, into a list of these terms .*)
1.2617 -fun term2list (t as (Const("HOL.divide",_) $ _ $ _)) = [t]
1.2618 - | term2list (t as (Const("Atools.pow",_) $ _ $ _)) =
1.2619 - [Const("HOL.divide",[HOLogic.realT,HOLogic.realT]--->HOLogic.realT) $
1.2620 - t $ Free("1",HOLogic.realT)
1.2621 - ]
1.2622 - | term2list (t as (Free(_,_))) =
1.2623 - [Const("HOL.divide",[HOLogic.realT,HOLogic.realT]--->HOLogic.realT) $
1.2624 - t $ Free("1",HOLogic.realT)
1.2625 - ]
1.2626 - | term2list (t as (Const("op *",_) $ _ $ _)) =
1.2627 - [Const("HOL.divide",[HOLogic.realT,HOLogic.realT]--->HOLogic.realT) $
1.2628 - t $ Free("1",HOLogic.realT)
1.2629 - ]
1.2630 - | term2list (Const("op +",_) $ t1 $ t2) = term2list(t1) @ term2list(t2)
1.2631 - | term2list (Const("op -",_) $ t1 $ t2) =
1.2632 - raise error ("RATIONALS_TERM2LIST_EXCEPTION: - not implemented yet")
1.2633 - | term2list _ = raise error ("RATIONALS_TERM2LIST_EXCEPTION: invalid term");
1.2634 -
1.2635 -(*.factors out the gcd of nominator and denominator:
1.2636 - a/b = (a' * gcd)/(b' * gcd), a,b,gcd are poly[2].*)
1.2637 -fun factout_p_ (thy:theory) t = SOME (step_cancel t,[]:term list);
1.2638 -fun factout_ (thy:theory) t = SOME (step_cancel_expanded t,[]:term list);
1.2639 -
1.2640 -(*.cancels a single fraction with normalform [2]
1.2641 - resulting in a canceled fraction [2], see factout_ .*)
1.2642 -fun cancel_p_ (thy:theory) t = (*WN.2.6.03 no rewrite -> NONE !*)
1.2643 - (let val (t',asm) = direct_cancel(*_expanded ... corrected MG.21.8.03*) t
1.2644 - in if t = t' then NONE else SOME (t',asm)
1.2645 - end) handle _ => NONE;
1.2646 -(*.the same as above with normalform [3]
1.2647 - val cancel_ :
1.2648 - theory -> (*10.02 unused *)
1.2649 - term -> (*fraction in normalform [3] *)
1.2650 - (term * (*fraction in normalform [3] *)
1.2651 - term list) (*casual asumptions in normalform [3] *)
1.2652 - option (*NONE: the function is not applicable *).*)
1.2653 -fun cancel_ (thy:theory) t = SOME (direct_cancel_expanded t) handle _ => NONE;
1.2654 -
1.2655 -(*.transforms sums of at least 2 fractions [3] to
1.2656 - sums with the least common multiple as nominator.*)
1.2657 -fun common_nominator_p_ (thy:theory) t =
1.2658 -((*writeln("### common_nominator_p_ called");*)
1.2659 - SOME (step_add_list_of_fractions(term2list(t))) handle _ => NONE
1.2660 -);
1.2661 -fun common_nominator_ (thy:theory) t =
1.2662 - SOME (step_add_list_of_fractions_exp(term2list(t))) handle _ => NONE;
1.2663 -
1.2664 -(*.add 2 or more fractions
1.2665 -val add_fraction_p_ :
1.2666 - theory -> (*10.02 unused *)
1.2667 - term -> (*2 or more fractions with normalform [2] *)
1.2668 - (term * (*one fraction with normalform [2] *)
1.2669 - term list) (*casual assumptions in normalform [2] WN0210???SK *)
1.2670 - option (*NONE: the function is not applicable *).*)
1.2671 -fun add_fraction_p_ (thy:theory) t =
1.2672 -(writeln("### add_fraction_p_ called");
1.2673 - (let val ts = term2list t
1.2674 - in if 1 < length ts
1.2675 - then SOME (add_list_of_fractions ts)
1.2676 - else NONE (*raise error ("RATIONALS_ADD_EXCEPTION: nothing to add")*)
1.2677 - end) handle _ => NONE
1.2678 -);
1.2679 -(*.same as add_fraction_p_ but with normalform [3].*)
1.2680 -(*SOME (step_add_list_of_fractions2(term2list(t))); *)
1.2681 -fun add_fraction_ (thy:theory) t =
1.2682 - if length(term2list(t))>1
1.2683 - then SOME (add_list_of_fractions_exp(term2list(t))) handle _ => NONE
1.2684 - else (*raise error ("RATIONALS_ADD_FRACTION_EXCEPTION: nothing to add")*)
1.2685 - NONE;
1.2686 -fun add_fraction_ (thy:theory) t =
1.2687 - (if 1 < length (term2list t)
1.2688 - then SOME (add_list_of_fractions_exp (term2list t))
1.2689 - else (*raise error ("RATIONALS_ADD_FRACTION_EXCEPTION: nothing to add")*)
1.2690 - NONE) handle _ => NONE;
1.2691 -
1.2692 -(*SOME (step_add_list_of_fractions2_exp(term2list(t))); *)
1.2693 -
1.2694 -(*. brings the term into a normal form .*)
1.2695 -fun norm_rational_ (thy:theory) t =
1.2696 - SOME (add_list_of_fractions(term2list(t))) handle _ => NONE;
1.2697 -fun norm_expanded_rat_ (thy:theory) t =
1.2698 - SOME (add_list_of_fractions_exp(term2list(t))) handle _ => NONE;
1.2699 -
1.2700 -
1.2701 -(*.evaluates conditions in calculate_Rational.*)
1.2702 -(*make local with FIXX@ME result:term *term list*)
1.2703 -val calc_rat_erls = prep_rls(
1.2704 - Rls {id = "calc_rat_erls", preconds = [], rew_ord = ("dummy_ord",dummy_ord),
1.2705 - erls = e_rls, srls = Erls, calc = [], (*asm_thm = [], *)
1.2706 - rules =
1.2707 - [Calc ("op =",eval_equal "#equal_"),
1.2708 - Calc ("Atools.is'_const",eval_const "#is_const_"),
1.2709 - Thm ("not_true",num_str not_true),
1.2710 - Thm ("not_false",num_str not_false)
1.2711 - ],
1.2712 - scr = EmptyScr});
1.2713 -
1.2714 -
1.2715 -(*.simplifies expressions with numerals;
1.2716 - does NOT rearrange the term by AC-rewriting; thus terms with variables
1.2717 - need to have constants to be commuted together respectively.*)
1.2718 -val calculate_Rational = prep_rls(
1.2719 - merge_rls "calculate_Rational"
1.2720 - (Rls {id = "divide", preconds = [], rew_ord = ("dummy_ord",dummy_ord),
1.2721 - erls = calc_rat_erls, srls = Erls, (*asm_thm = [],*)
1.2722 - calc = [],
1.2723 - rules =
1.2724 - [Calc ("HOL.divide" ,eval_cancel "#divide_"),
1.2725 -
1.2726 - Thm ("sym_real_minus_divide_eq",
1.2727 - num_str (real_minus_divide_eq RS sym)),
1.2728 - (*SYM - ?x / ?y = - (?x / ?y) may come from subst*)
1.2729 -
1.2730 - Thm ("rat_add",num_str rat_add),
1.2731 - (*"[| a is_const; b is_const; c is_const; d is_const |] ==> \
1.2732 - \"a / c + b / d = (a * d) / (c * d) + (b * c ) / (d * c)"*)
1.2733 - Thm ("rat_add1",num_str rat_add1),
1.2734 - (*"[| a is_const; b is_const; c is_const |] ==> \
1.2735 - \"a / c + b / c = (a + b) / c"*)
1.2736 - Thm ("rat_add2",num_str rat_add2),
1.2737 - (*"[| ?a is_const; ?b is_const; ?c is_const |] ==> \
1.2738 - \?a / ?c + ?b = (?a + ?b * ?c) / ?c"*)
1.2739 - Thm ("rat_add3",num_str rat_add3),
1.2740 - (*"[| a is_const; b is_const; c is_const |] ==> \
1.2741 - \"a + b / c = (a * c) / c + b / c"\
1.2742 - \.... is_const to be omitted here FIXME*)
1.2743 -
1.2744 - Thm ("rat_mult",num_str rat_mult),
1.2745 - (*a / b * (c / d) = a * c / (b * d)*)
1.2746 - Thm ("real_times_divide1_eq",num_str real_times_divide1_eq),
1.2747 - (*?x * (?y / ?z) = ?x * ?y / ?z*)
1.2748 - Thm ("real_times_divide2_eq",num_str real_times_divide2_eq),
1.2749 - (*?y / ?z * ?x = ?y * ?x / ?z*)
1.2750 -
1.2751 - Thm ("real_divide_divide1",num_str real_divide_divide1),
1.2752 - (*"?y ~= 0 ==> ?u / ?v / (?y / ?z) = ?u / ?v * (?z / ?y)"*)
1.2753 - Thm ("real_divide_divide2_eq",num_str real_divide_divide2_eq),
1.2754 - (*"?x / ?y / ?z = ?x / (?y * ?z)"*)
1.2755 -
1.2756 - Thm ("rat_power", num_str rat_power),
1.2757 - (*"(?a / ?b) ^^^ ?n = ?a ^^^ ?n / ?b ^^^ ?n"*)
1.2758 -
1.2759 - Thm ("mult_cross",num_str mult_cross),
1.2760 - (*"[| b ~= 0; d ~= 0 |] ==> (a / b = c / d) = (a * d = b * c)*)
1.2761 - Thm ("mult_cross1",num_str mult_cross1),
1.2762 - (*" b ~= 0 ==> (a / b = c ) = (a = b * c)*)
1.2763 - Thm ("mult_cross2",num_str mult_cross2)
1.2764 - (*" d ~= 0 ==> (a = c / d) = (a * d = c)*)
1.2765 - ], scr = EmptyScr})
1.2766 - calculate_Poly);
1.2767 -
1.2768 -
1.2769 -(*("is_expanded", ("Rational.is'_expanded", eval_is_expanded ""))*)
1.2770 -fun eval_is_expanded (thmid:string) _
1.2771 - (t as (Const("Rational.is'_expanded", _) $ arg)) thy =
1.2772 - if is_expanded arg
1.2773 - then SOME (mk_thmid thmid ""
1.2774 - ((Syntax.string_of_term (thy2ctxt thy)) arg) "",
1.2775 - Trueprop $ (mk_equality (t, HOLogic.true_const)))
1.2776 - else SOME (mk_thmid thmid ""
1.2777 - ((Syntax.string_of_term (thy2ctxt thy)) arg) "",
1.2778 - Trueprop $ (mk_equality (t, HOLogic.false_const)))
1.2779 - | eval_is_expanded _ _ _ _ = NONE;
1.2780 -
1.2781 -val rational_erls =
1.2782 - merge_rls "rational_erls" calculate_Rational
1.2783 - (append_rls "is_expanded" Atools_erls
1.2784 - [Calc ("Rational.is'_expanded", eval_is_expanded "")
1.2785 - ]);
1.2786 -
1.2787 -
1.2788 -
1.2789 -(*.3 'reverse-rewrite-sets' for symbolic computation on rationals:
1.2790 - =================================================================
1.2791 - A[2] 'cancel_p': .
1.2792 - A[3] 'cancel': .
1.2793 - B[2] 'common_nominator_p': transforms summands in a term [2]
1.2794 - to fractions with the (least) common multiple as nominator.
1.2795 - B[3] 'norm_rational': normalizes arbitrary algebraic terms (without
1.2796 - radicals and transzendental functions) to one canceled fraction,
1.2797 - nominator and denominator in polynomial form.
1.2798 -
1.2799 -In order to meet isac's requirements for interactive and stepwise calculation,
1.2800 -each 'reverse-rewerite-set' consists of an initialization for the interpreter
1.2801 -state and of 4 functions, each of which employs rewriting as much as possible.
1.2802 -The signature of these functions are the same in each 'reverse-rewrite-set'
1.2803 -respectively.*)
1.2804 -
1.2805 -(* ************************************************************************* *)
1.2806 -
1.2807 -
1.2808 -local(*. cancel_p
1.2809 -------------------------
1.2810 -cancels a single fraction consisting of two (uni- or multivariate)
1.2811 -polynomials WN0609???SK[2] into another such a fraction; examples:
1.2812 -
1.2813 - a^2 + -1*b^2 a + b
1.2814 - -------------------- = ---------
1.2815 - a^2 + -2*a*b + b^2 a + -1*b
1.2816 -
1.2817 - a^2 a
1.2818 - --- = ---
1.2819 - a 1
1.2820 -
1.2821 -Remark: the reverse ruleset does _NOT_ work properly with other input !.*)
1.2822 -(*WN020824 wir werden "uberlegen, wie wir ungeeignete inputs zur"uckweisen*)
1.2823 -
1.2824 -val {rules, rew_ord=(_,ro),...} =
1.2825 - rep_rls (assoc_rls "make_polynomial");
1.2826 -(*WN060829 ... make_deriv does not terminate with 1st expl above,
1.2827 - see rational.sml --- investigate rulesets for cancel_p ---*)
1.2828 -val {rules, rew_ord=(_,ro),...} =
1.2829 - rep_rls (assoc_rls "rev_rew_p");
1.2830 -
1.2831 -val thy = Rational.thy;
1.2832 -
1.2833 -(*.init_state = fn : term -> istate
1.2834 -initialzies the state of the script interpreter. The state is:
1.2835 -
1.2836 -type rrlsstate = (*state for reverse rewriting*)
1.2837 - (term * (*the current formula*)
1.2838 - term * (*the final term*)
1.2839 - rule list (*'reverse rule list' (#)*)
1.2840 - list * (*may be serveral, eg. in norm_rational*)
1.2841 - (rule * (*Thm (+ Thm generated from Calc) resulting in ...*)
1.2842 - (term * (*... rewrite with ...*)
1.2843 - term list)) (*... assumptions*)
1.2844 - list); (*derivation from given term to normalform
1.2845 - in reverse order with sym_thm;
1.2846 - (#) could be extracted from here by (map #1)*).*)
1.2847 -(* val {rules, rew_ord=(_,ro),...} =
1.2848 - rep_rls (assoc_rls "rev_rew_p") (*USE ALWAYS, SEE val cancel_p*);
1.2849 - val (thy, eval_rls, ro) =(Rational.thy, Atools_erls, ro) (*..val cancel_p*);
1.2850 - val t = t;
1.2851 - *)
1.2852 -fun init_state thy eval_rls ro t =
1.2853 - let val SOME (t',_) = factout_p_ thy t
1.2854 - val SOME (t'',asm) = cancel_p_ thy t
1.2855 - val der = reverse_deriv thy eval_rls rules ro NONE t'
1.2856 - val der = der @ [(Thm ("real_mult_div_cancel2",
1.2857 - num_str real_mult_div_cancel2),
1.2858 - (t'',asm))]
1.2859 - val rs = (distinct_Thm o (map #1)) der
1.2860 - val rs = filter_out (eq_Thms ["sym_real_add_zero_left",
1.2861 - "sym_real_mult_0",
1.2862 - "sym_real_mult_1"
1.2863 - (*..insufficient,eg.make_Polynomial*)])rs
1.2864 - in (t,t'',[rs(*here only _ONE_ to ease locate_rule*)],der) end;
1.2865 -
1.2866 -(*.locate_rule = fn : rule list -> term -> rule
1.2867 - -> (rule * (term * term list) option) list.
1.2868 - checks a rule R for being a cancel-rule, and if it is,
1.2869 - then return the list of rules (+ the terms they are rewriting to)
1.2870 - which need to be applied before R should be applied.
1.2871 - precondition: the rule is applicable to the argument-term.
1.2872 -arguments:
1.2873 - rule list: the reverse rule list
1.2874 - -> term : ... to which the rule shall be applied
1.2875 - -> rule : ... to be applied to term
1.2876 -value:
1.2877 - -> (rule : a rule rewriting to ...
1.2878 - * (term : ... the resulting term ...
1.2879 - * term list): ... with the assumptions ( //#0).
1.2880 - ) list : there may be several such rules;
1.2881 - the list is empty, if the rule has nothing to do
1.2882 - with cancelation.*)
1.2883 -(* val () = ();
1.2884 - *)
1.2885 -fun locate_rule thy eval_rls ro [rs] t r =
1.2886 - if (id_of_thm r) mem (map (id_of_thm)) rs
1.2887 - then let val ropt =
1.2888 - rewrite_ thy ro eval_rls true (thm_of_thm r) t;
1.2889 - in case ropt of
1.2890 - SOME ta => [(r, ta)]
1.2891 - | NONE => (writeln("### locate_rule: rewrite "^
1.2892 - (id_of_thm r)^" "^(term2str t)^" = NONE");
1.2893 - []) end
1.2894 - else (writeln("### locate_rule: "^(id_of_thm r)^" not mem rrls");[])
1.2895 - | locate_rule _ _ _ _ _ _ =
1.2896 - raise error ("locate_rule: doesnt match rev-sets in istate");
1.2897 -
1.2898 -(*.next_rule = fn : rule list -> term -> rule option
1.2899 - for a given term return the next rules to be done for cancelling.
1.2900 -arguments:
1.2901 - rule list : the reverse rule list
1.2902 - term : the term for which ...
1.2903 -value:
1.2904 - -> rule option: ... this rule is appropriate for cancellation;
1.2905 - there may be no such rule (if the term is canceled already.*)
1.2906 -(* val thy = Rational.thy;
1.2907 - val Rrls {rew_ord=(_,ro),...} = cancel;
1.2908 - val ([rs],t) = (rss,f);
1.2909 - next_rule thy eval_rls ro [rs] t;(*eval fun next_rule ... before!*)
1.2910 -
1.2911 - val (thy, [rs]) = (Rational.thy, revsets);
1.2912 - val Rrls {rew_ord=(_,ro),...} = cancel;
1.2913 - nex [rs] t;
1.2914 - *)
1.2915 -fun next_rule thy eval_rls ro [rs] t =
1.2916 - let val der = make_deriv thy eval_rls rs ro NONE t;
1.2917 - in case der of
1.2918 -(* val (_,r,_)::_ = der;
1.2919 - *)
1.2920 - (_,r,_)::_ => SOME r
1.2921 - | _ => NONE
1.2922 - end
1.2923 - | next_rule _ _ _ _ _ =
1.2924 - raise error ("next_rule: doesnt match rev-sets in istate");
1.2925 -
1.2926 -(*.val attach_form = f : rule list -> term -> term
1.2927 - -> (rule * (term * term list)) list
1.2928 - checks an input term TI, if it may belong to a current cancellation, by
1.2929 - trying to derive it from the given term TG.
1.2930 -arguments:
1.2931 - term : TG, the last one in the cancellation agreed upon by user + math-eng
1.2932 - -> term: TI, the next one input by the user
1.2933 -value:
1.2934 - -> (rule : the rule to be applied in order to reach TI
1.2935 - * (term : ... obtained by applying the rule ...
1.2936 - * term list): ... and the respective assumptions.
1.2937 - ) list : there may be several such rules;
1.2938 - the list is empty, if the users term does not belong
1.2939 - to a cancellation of the term last agreed upon.*)
1.2940 -fun attach_form (_:rule list list) (_:term) (_:term) = (*still missing*)
1.2941 - []:(rule * (term * term list)) list;
1.2942 -
1.2943 -in
1.2944 -
1.2945 -val cancel_p =
1.2946 - Rrls {id = "cancel_p", prepat=[],
1.2947 - rew_ord=("ord_make_polynomial",
1.2948 - ord_make_polynomial false Rational.thy),
1.2949 - erls = rational_erls,
1.2950 - calc = [("PLUS" ,("op +" ,eval_binop "#add_")),
1.2951 - ("TIMES" ,("op *" ,eval_binop "#mult_")),
1.2952 - ("DIVIDE" ,("HOL.divide" ,eval_cancel "#divide_")),
1.2953 - ("POWER" ,("Atools.pow" ,eval_binop "#power_"))],
1.2954 - (*asm_thm=[("real_mult_div_cancel2","")],*)
1.2955 - scr=Rfuns {init_state = init_state thy Atools_erls ro,
1.2956 - normal_form = cancel_p_ thy,
1.2957 - locate_rule = locate_rule thy Atools_erls ro,
1.2958 - next_rule = next_rule thy Atools_erls ro,
1.2959 - attach_form = attach_form}}
1.2960 -end;(*local*)
1.2961 -
1.2962 -
1.2963 -local(*.ad (1) 'cancel'
1.2964 -------------------------------
1.2965 -cancels a single fraction consisting of two (uni- or multivariate)
1.2966 -polynomials WN0609???SK[3] into another such a fraction; examples:
1.2967 -
1.2968 - a^2 - b^2 a + b
1.2969 - -------------------- = ---------
1.2970 - a^2 - 2*a*b + b^2 a - *b
1.2971 -
1.2972 -Remark: the reverse ruleset does _NOT_ work properly with other input !.*)
1.2973 -(*WN 24.8.02: wir werden "uberlegen, wie wir ungeeignete inputs zur"uckweisen*)
1.2974 -
1.2975 -(*
1.2976 -val SOME (Rls {rules=rules,rew_ord=(_,ro),...}) =
1.2977 - assoc'(!ruleset',"expand_binoms");
1.2978 -*)
1.2979 -val {rules=rules,rew_ord=(_,ro),...} =
1.2980 - rep_rls (assoc_rls "expand_binoms");
1.2981 -val thy = Rational.thy;
1.2982 -
1.2983 -fun init_state thy eval_rls ro t =
1.2984 - let val SOME (t',_) = factout_ thy t;
1.2985 - val SOME (t'',asm) = cancel_ thy t;
1.2986 - val der = reverse_deriv thy eval_rls rules ro NONE t';
1.2987 - val der = der @ [(Thm ("real_mult_div_cancel2",
1.2988 - num_str real_mult_div_cancel2),
1.2989 - (t'',asm))]
1.2990 - val rs = map #1 der;
1.2991 - in (t,t'',[rs],der) end;
1.2992 -
1.2993 -fun locate_rule thy eval_rls ro [rs] t r =
1.2994 - if (id_of_thm r) mem (map (id_of_thm)) rs
1.2995 - then let val ropt =
1.2996 - rewrite_ thy ro eval_rls true (thm_of_thm r) t;
1.2997 - in case ropt of
1.2998 - SOME ta => [(r, ta)]
1.2999 - | NONE => (writeln("### locate_rule: rewrite "^
1.3000 - (id_of_thm r)^" "^(term2str t)^" = NONE");
1.3001 - []) end
1.3002 - else (writeln("### locate_rule: "^(id_of_thm r)^" not mem rrls");[])
1.3003 - | locate_rule _ _ _ _ _ _ =
1.3004 - raise error ("locate_rule: doesnt match rev-sets in istate");
1.3005 -
1.3006 -fun next_rule thy eval_rls ro [rs] t =
1.3007 - let val der = make_deriv thy eval_rls rs ro NONE t;
1.3008 - in case der of
1.3009 -(* val (_,r,_)::_ = der;
1.3010 - *)
1.3011 - (_,r,_)::_ => SOME r
1.3012 - | _ => NONE
1.3013 - end
1.3014 - | next_rule _ _ _ _ _ =
1.3015 - raise error ("next_rule: doesnt match rev-sets in istate");
1.3016 -
1.3017 -fun attach_form (_:rule list list) (_:term) (_:term) = (*still missing*)
1.3018 - []:(rule * (term * term list)) list;
1.3019 -
1.3020 -val pat = (term_of o the o (parse thy)) "?r/?s";
1.3021 -val pre1 = (term_of o the o (parse thy)) "?r is_expanded";
1.3022 -val pre2 = (term_of o the o (parse thy)) "?s is_expanded";
1.3023 -val prepat = [([pre1, pre2], pat)];
1.3024 -
1.3025 -in
1.3026 -
1.3027 -
1.3028 -val cancel =
1.3029 - Rrls {id = "cancel", prepat=prepat,
1.3030 - rew_ord=("ord_make_polynomial",
1.3031 - ord_make_polynomial false Rational.thy),
1.3032 - erls = rational_erls,
1.3033 - calc = [("PLUS" ,("op +" ,eval_binop "#add_")),
1.3034 - ("TIMES" ,("op *" ,eval_binop "#mult_")),
1.3035 - ("DIVIDE" ,("HOL.divide" ,eval_cancel "#divide_")),
1.3036 - ("POWER" ,("Atools.pow" ,eval_binop "#power_"))],
1.3037 - scr=Rfuns {init_state = init_state thy Atools_erls ro,
1.3038 - normal_form = cancel_ thy,
1.3039 - locate_rule = locate_rule thy Atools_erls ro,
1.3040 - next_rule = next_rule thy Atools_erls ro,
1.3041 - attach_form = attach_form}}
1.3042 -end;(*local*)
1.3043 -
1.3044 -
1.3045 -
1.3046 -local(*.ad [2] 'common_nominator_p'
1.3047 ----------------------------------
1.3048 -FIXME Beschreibung .*)
1.3049 -
1.3050 -
1.3051 -val {rules=rules,rew_ord=(_,ro),...} =
1.3052 - rep_rls (assoc_rls "make_polynomial");
1.3053 -(*WN060829 ... make_deriv does not terminate with 1st expl above,
1.3054 - see rational.sml --- investigate rulesets for cancel_p ---*)
1.3055 -val {rules, rew_ord=(_,ro),...} =
1.3056 - rep_rls (assoc_rls "rev_rew_p");
1.3057 -val thy = Rational.thy;
1.3058 -
1.3059 -
1.3060 -(*.common_nominator_p_ = fn : theory -> term -> (term * term list) option
1.3061 - as defined above*)
1.3062 -
1.3063 -(*.init_state = fn : term -> istate
1.3064 -initialzies the state of the interactive interpreter. The state is:
1.3065 -
1.3066 -type rrlsstate = (*state for reverse rewriting*)
1.3067 - (term * (*the current formula*)
1.3068 - term * (*the final term*)
1.3069 - rule list (*'reverse rule list' (#)*)
1.3070 - list * (*may be serveral, eg. in norm_rational*)
1.3071 - (rule * (*Thm (+ Thm generated from Calc) resulting in ...*)
1.3072 - (term * (*... rewrite with ...*)
1.3073 - term list)) (*... assumptions*)
1.3074 - list); (*derivation from given term to normalform
1.3075 - in reverse order with sym_thm;
1.3076 - (#) could be extracted from here by (map #1)*).*)
1.3077 -fun init_state thy eval_rls ro t =
1.3078 - let val SOME (t',_) = common_nominator_p_ thy t;
1.3079 - val SOME (t'',asm) = add_fraction_p_ thy t;
1.3080 - val der = reverse_deriv thy eval_rls rules ro NONE t';
1.3081 - val der = der @ [(Thm ("real_mult_div_cancel2",
1.3082 - num_str real_mult_div_cancel2),
1.3083 - (t'',asm))]
1.3084 - val rs = (distinct_Thm o (map #1)) der;
1.3085 - val rs = filter_out (eq_Thms ["sym_real_add_zero_left",
1.3086 - "sym_real_mult_0",
1.3087 - "sym_real_mult_1"]) rs;
1.3088 - in (t,t'',[rs(*here only _ONE_*)],der) end;
1.3089 -
1.3090 -(* use"knowledge/Rational.ML";
1.3091 - *)
1.3092 -
1.3093 -(*.locate_rule = fn : rule list -> term -> rule
1.3094 - -> (rule * (term * term list) option) list.
1.3095 - checks a rule R for being a cancel-rule, and if it is,
1.3096 - then return the list of rules (+ the terms they are rewriting to)
1.3097 - which need to be applied before R should be applied.
1.3098 - precondition: the rule is applicable to the argument-term.
1.3099 -arguments:
1.3100 - rule list: the reverse rule list
1.3101 - -> term : ... to which the rule shall be applied
1.3102 - -> rule : ... to be applied to term
1.3103 -value:
1.3104 - -> (rule : a rule rewriting to ...
1.3105 - * (term : ... the resulting term ...
1.3106 - * term list): ... with the assumptions ( //#0).
1.3107 - ) list : there may be several such rules;
1.3108 - the list is empty, if the rule has nothing to do
1.3109 - with cancelation.*)
1.3110 -(* val () = ();
1.3111 - *)
1.3112 -fun locate_rule thy eval_rls ro [rs] t r =
1.3113 - if (id_of_thm r) mem (map (id_of_thm)) rs
1.3114 - then let val ropt =
1.3115 - rewrite_ thy ro eval_rls true (thm_of_thm r) t;
1.3116 - in case ropt of
1.3117 - SOME ta => [(r, ta)]
1.3118 - | NONE => (writeln("### locate_rule: rewrite "^
1.3119 - (id_of_thm r)^" "^(term2str t)^" = NONE");
1.3120 - []) end
1.3121 - else (writeln("### locate_rule: "^(id_of_thm r)^" not mem rrls");[])
1.3122 - | locate_rule _ _ _ _ _ _ =
1.3123 - raise error ("locate_rule: doesnt match rev-sets in istate");
1.3124 -
1.3125 -(*.next_rule = fn : rule list -> term -> rule option
1.3126 - for a given term return the next rules to be done for cancelling.
1.3127 -arguments:
1.3128 - rule list : the reverse rule list
1.3129 - term : the term for which ...
1.3130 -value:
1.3131 - -> rule option: ... this rule is appropriate for cancellation;
1.3132 - there may be no such rule (if the term is canceled already.*)
1.3133 -(* val thy = Rational.thy;
1.3134 - val Rrls {rew_ord=(_,ro),...} = cancel;
1.3135 - val ([rs],t) = (rss,f);
1.3136 - next_rule thy eval_rls ro [rs] t;(*eval fun next_rule ... before!*)
1.3137 -
1.3138 - val (thy, [rs]) = (Rational.thy, revsets);
1.3139 - val Rrls {rew_ord=(_,ro),...} = cancel;
1.3140 - nex [rs] t;
1.3141 - *)
1.3142 -fun next_rule thy eval_rls ro [rs] t =
1.3143 - let val der = make_deriv thy eval_rls rs ro NONE t;
1.3144 - in case der of
1.3145 -(* val (_,r,_)::_ = der;
1.3146 - *)
1.3147 - (_,r,_)::_ => SOME r
1.3148 - | _ => NONE
1.3149 - end
1.3150 - | next_rule _ _ _ _ _ =
1.3151 - raise error ("next_rule: doesnt match rev-sets in istate");
1.3152 -
1.3153 -(*.val attach_form = f : rule list -> term -> term
1.3154 - -> (rule * (term * term list)) list
1.3155 - checks an input term TI, if it may belong to a current cancellation, by
1.3156 - trying to derive it from the given term TG.
1.3157 -arguments:
1.3158 - term : TG, the last one in the cancellation agreed upon by user + math-eng
1.3159 - -> term: TI, the next one input by the user
1.3160 -value:
1.3161 - -> (rule : the rule to be applied in order to reach TI
1.3162 - * (term : ... obtained by applying the rule ...
1.3163 - * term list): ... and the respective assumptions.
1.3164 - ) list : there may be several such rules;
1.3165 - the list is empty, if the users term does not belong
1.3166 - to a cancellation of the term last agreed upon.*)
1.3167 -fun attach_form (_:rule list list) (_:term) (_:term) = (*still missing*)
1.3168 - []:(rule * (term * term list)) list;
1.3169 -
1.3170 -val pat0 = (term_of o the o (parse thy)) "?r/?s+?u/?v";
1.3171 -val pat1 = (term_of o the o (parse thy)) "?r/?s+?u ";
1.3172 -val pat2 = (term_of o the o (parse thy)) "?r +?u/?v";
1.3173 -val prepat = [([HOLogic.true_const], pat0),
1.3174 - ([HOLogic.true_const], pat1),
1.3175 - ([HOLogic.true_const], pat2)];
1.3176 -
1.3177 -in
1.3178 -
1.3179 -(*11.02 schnelle L"osung f"ur RL: Bruch auch gek"urzt;
1.3180 - besser w"are: auf 1 gemeinsamen Bruchstrich, Nenner und Z"ahler unvereinfacht
1.3181 - dh. wie common_nominator_p_, aber auf 1 Bruchstrich*)
1.3182 -val common_nominator_p =
1.3183 - Rrls {id = "common_nominator_p", prepat=prepat,
1.3184 - rew_ord=("ord_make_polynomial",
1.3185 - ord_make_polynomial false Rational.thy),
1.3186 - erls = rational_erls,
1.3187 - calc = [("PLUS" ,("op +" ,eval_binop "#add_")),
1.3188 - ("TIMES" ,("op *" ,eval_binop "#mult_")),
1.3189 - ("DIVIDE" ,("HOL.divide" ,eval_cancel "#divide_")),
1.3190 - ("POWER" ,("Atools.pow" ,eval_binop "#power_"))],
1.3191 - scr=Rfuns {init_state = init_state thy Atools_erls ro,
1.3192 - normal_form = add_fraction_p_ thy,(*FIXME.WN0211*)
1.3193 - locate_rule = locate_rule thy Atools_erls ro,
1.3194 - next_rule = next_rule thy Atools_erls ro,
1.3195 - attach_form = attach_form}}
1.3196 -end;(*local*)
1.3197 -
1.3198 -
1.3199 -local(*.ad [2] 'common_nominator'
1.3200 ----------------------------------
1.3201 -FIXME Beschreibung .*)
1.3202 -
1.3203 -
1.3204 -val {rules=rules,rew_ord=(_,ro),...} =
1.3205 - rep_rls (assoc_rls "make_polynomial");
1.3206 -val thy = Rational.thy;
1.3207 -
1.3208 -
1.3209 -(*.common_nominator_ = fn : theory -> term -> (term * term list) option
1.3210 - as defined above*)
1.3211 -
1.3212 -(*.init_state = fn : term -> istate
1.3213 -initialzies the state of the interactive interpreter. The state is:
1.3214 -
1.3215 -type rrlsstate = (*state for reverse rewriting*)
1.3216 - (term * (*the current formula*)
1.3217 - term * (*the final term*)
1.3218 - rule list (*'reverse rule list' (#)*)
1.3219 - list * (*may be serveral, eg. in norm_rational*)
1.3220 - (rule * (*Thm (+ Thm generated from Calc) resulting in ...*)
1.3221 - (term * (*... rewrite with ...*)
1.3222 - term list)) (*... assumptions*)
1.3223 - list); (*derivation from given term to normalform
1.3224 - in reverse order with sym_thm;
1.3225 - (#) could be extracted from here by (map #1)*).*)
1.3226 -fun init_state thy eval_rls ro t =
1.3227 - let val SOME (t',_) = common_nominator_ thy t;
1.3228 - val SOME (t'',asm) = add_fraction_ thy t;
1.3229 - val der = reverse_deriv thy eval_rls rules ro NONE t';
1.3230 - val der = der @ [(Thm ("real_mult_div_cancel2",
1.3231 - num_str real_mult_div_cancel2),
1.3232 - (t'',asm))]
1.3233 - val rs = (distinct_Thm o (map #1)) der;
1.3234 - val rs = filter_out (eq_Thms ["sym_real_add_zero_left",
1.3235 - "sym_real_mult_0",
1.3236 - "sym_real_mult_1"]) rs;
1.3237 - in (t,t'',[rs(*here only _ONE_*)],der) end;
1.3238 -
1.3239 -(* use"knowledge/Rational.ML";
1.3240 - *)
1.3241 -
1.3242 -(*.locate_rule = fn : rule list -> term -> rule
1.3243 - -> (rule * (term * term list) option) list.
1.3244 - checks a rule R for being a cancel-rule, and if it is,
1.3245 - then return the list of rules (+ the terms they are rewriting to)
1.3246 - which need to be applied before R should be applied.
1.3247 - precondition: the rule is applicable to the argument-term.
1.3248 -arguments:
1.3249 - rule list: the reverse rule list
1.3250 - -> term : ... to which the rule shall be applied
1.3251 - -> rule : ... to be applied to term
1.3252 -value:
1.3253 - -> (rule : a rule rewriting to ...
1.3254 - * (term : ... the resulting term ...
1.3255 - * term list): ... with the assumptions ( //#0).
1.3256 - ) list : there may be several such rules;
1.3257 - the list is empty, if the rule has nothing to do
1.3258 - with cancelation.*)
1.3259 -(* val () = ();
1.3260 - *)
1.3261 -fun locate_rule thy eval_rls ro [rs] t r =
1.3262 - if (id_of_thm r) mem (map (id_of_thm)) rs
1.3263 - then let val ropt =
1.3264 - rewrite_ thy ro eval_rls true (thm_of_thm r) t;
1.3265 - in case ropt of
1.3266 - SOME ta => [(r, ta)]
1.3267 - | NONE => (writeln("### locate_rule: rewrite "^
1.3268 - (id_of_thm r)^" "^(term2str t)^" = NONE");
1.3269 - []) end
1.3270 - else (writeln("### locate_rule: "^(id_of_thm r)^" not mem rrls");[])
1.3271 - | locate_rule _ _ _ _ _ _ =
1.3272 - raise error ("locate_rule: doesnt match rev-sets in istate");
1.3273 -
1.3274 -(*.next_rule = fn : rule list -> term -> rule option
1.3275 - for a given term return the next rules to be done for cancelling.
1.3276 -arguments:
1.3277 - rule list : the reverse rule list
1.3278 - term : the term for which ...
1.3279 -value:
1.3280 - -> rule option: ... this rule is appropriate for cancellation;
1.3281 - there may be no such rule (if the term is canceled already.*)
1.3282 -(* val thy = Rational.thy;
1.3283 - val Rrls {rew_ord=(_,ro),...} = cancel;
1.3284 - val ([rs],t) = (rss,f);
1.3285 - next_rule thy eval_rls ro [rs] t;(*eval fun next_rule ... before!*)
1.3286 -
1.3287 - val (thy, [rs]) = (Rational.thy, revsets);
1.3288 - val Rrls {rew_ord=(_,ro),...} = cancel_p;
1.3289 - nex [rs] t;
1.3290 - *)
1.3291 -fun next_rule thy eval_rls ro [rs] t =
1.3292 - let val der = make_deriv thy eval_rls rs ro NONE t;
1.3293 - in case der of
1.3294 -(* val (_,r,_)::_ = der;
1.3295 - *)
1.3296 - (_,r,_)::_ => SOME r
1.3297 - | _ => NONE
1.3298 - end
1.3299 - | next_rule _ _ _ _ _ =
1.3300 - raise error ("next_rule: doesnt match rev-sets in istate");
1.3301 -
1.3302 -(*.val attach_form = f : rule list -> term -> term
1.3303 - -> (rule * (term * term list)) list
1.3304 - checks an input term TI, if it may belong to a current cancellation, by
1.3305 - trying to derive it from the given term TG.
1.3306 -arguments:
1.3307 - term : TG, the last one in the cancellation agreed upon by user + math-eng
1.3308 - -> term: TI, the next one input by the user
1.3309 -value:
1.3310 - -> (rule : the rule to be applied in order to reach TI
1.3311 - * (term : ... obtained by applying the rule ...
1.3312 - * term list): ... and the respective assumptions.
1.3313 - ) list : there may be several such rules;
1.3314 - the list is empty, if the users term does not belong
1.3315 - to a cancellation of the term last agreed upon.*)
1.3316 -fun attach_form (_:rule list list) (_:term) (_:term) = (*still missing*)
1.3317 - []:(rule * (term * term list)) list;
1.3318 -
1.3319 -val pat0 = (term_of o the o (parse thy)) "?r/?s+?u/?v";
1.3320 -val pat01 = (term_of o the o (parse thy)) "?r/?s-?u/?v";
1.3321 -val pat1 = (term_of o the o (parse thy)) "?r/?s+?u ";
1.3322 -val pat11 = (term_of o the o (parse thy)) "?r/?s-?u ";
1.3323 -val pat2 = (term_of o the o (parse thy)) "?r +?u/?v";
1.3324 -val pat21 = (term_of o the o (parse thy)) "?r -?u/?v";
1.3325 -val prepat = [([HOLogic.true_const], pat0),
1.3326 - ([HOLogic.true_const], pat01),
1.3327 - ([HOLogic.true_const], pat1),
1.3328 - ([HOLogic.true_const], pat11),
1.3329 - ([HOLogic.true_const], pat2),
1.3330 - ([HOLogic.true_const], pat21)];
1.3331 -
1.3332 -
1.3333 -in
1.3334 -
1.3335 -val common_nominator =
1.3336 - Rrls {id = "common_nominator", prepat=prepat,
1.3337 - rew_ord=("ord_make_polynomial",
1.3338 - ord_make_polynomial false Rational.thy),
1.3339 - erls = rational_erls,
1.3340 - calc = [("PLUS" ,("op +" ,eval_binop "#add_")),
1.3341 - ("TIMES" ,("op *" ,eval_binop "#mult_")),
1.3342 - ("DIVIDE" ,("HOL.divide" ,eval_cancel "#divide_")),
1.3343 - ("POWER" ,("Atools.pow" ,eval_binop "#power_"))],
1.3344 - (*asm_thm=[("real_mult_div_cancel2","")],*)
1.3345 - scr=Rfuns {init_state = init_state thy Atools_erls ro,
1.3346 - normal_form = add_fraction_ (*NOT common_nominator_*) thy,
1.3347 - locate_rule = locate_rule thy Atools_erls ro,
1.3348 - next_rule = next_rule thy Atools_erls ro,
1.3349 - attach_form = attach_form}}
1.3350 -
1.3351 -end;(*local*)
1.3352 -
1.3353 -
1.3354 -(*##*)
1.3355 -end;(*struct*)
1.3356 -
1.3357 -open RationalI;
1.3358 -(*##*)
1.3359 -
1.3360 -(*.the expression contains + - * ^ / only ?.*)
1.3361 -fun is_ratpolyexp (Free _) = true
1.3362 - | is_ratpolyexp (Const ("op +",_) $ Free _ $ Free _) = true
1.3363 - | is_ratpolyexp (Const ("op -",_) $ Free _ $ Free _) = true
1.3364 - | is_ratpolyexp (Const ("op *",_) $ Free _ $ Free _) = true
1.3365 - | is_ratpolyexp (Const ("Atools.pow",_) $ Free _ $ Free _) = true
1.3366 - | is_ratpolyexp (Const ("HOL.divide",_) $ Free _ $ Free _) = true
1.3367 - | is_ratpolyexp (Const ("op +",_) $ t1 $ t2) =
1.3368 - ((is_ratpolyexp t1) andalso (is_ratpolyexp t2))
1.3369 - | is_ratpolyexp (Const ("op -",_) $ t1 $ t2) =
1.3370 - ((is_ratpolyexp t1) andalso (is_ratpolyexp t2))
1.3371 - | is_ratpolyexp (Const ("op *",_) $ t1 $ t2) =
1.3372 - ((is_ratpolyexp t1) andalso (is_ratpolyexp t2))
1.3373 - | is_ratpolyexp (Const ("Atools.pow",_) $ t1 $ t2) =
1.3374 - ((is_ratpolyexp t1) andalso (is_ratpolyexp t2))
1.3375 - | is_ratpolyexp (Const ("HOL.divide",_) $ t1 $ t2) =
1.3376 - ((is_ratpolyexp t1) andalso (is_ratpolyexp t2))
1.3377 - | is_ratpolyexp _ = false;
1.3378 -
1.3379 -(*("is_ratpolyexp", ("Rational.is'_ratpolyexp", eval_is_ratpolyexp ""))*)
1.3380 -fun eval_is_ratpolyexp (thmid:string) _
1.3381 - (t as (Const("Rational.is'_ratpolyexp", _) $ arg)) thy =
1.3382 - if is_ratpolyexp arg
1.3383 - then SOME (mk_thmid thmid ""
1.3384 - ((Syntax.string_of_term (thy2ctxt thy)) arg) "",
1.3385 - Trueprop $ (mk_equality (t, HOLogic.true_const)))
1.3386 - else SOME (mk_thmid thmid ""
1.3387 - ((Syntax.string_of_term (thy2ctxt thy)) arg) "",
1.3388 - Trueprop $ (mk_equality (t, HOLogic.false_const)))
1.3389 - | eval_is_ratpolyexp _ _ _ _ = NONE;
1.3390 -
1.3391 -
1.3392 -
1.3393 -(*-------------------18.3.03 --> struct <-----------vvv--*)
1.3394 -val add_fractions_p = common_nominator_p; (*FIXXXME:eilig f"ur norm_Rational*)
1.3395 -
1.3396 -(*.discard binary minus, shift unary minus into -1*;
1.3397 - unary minus before numerals are put into the numeral by parsing;
1.3398 - contains absolute minimum of thms for context in norm_Rational .*)
1.3399 -val discard_minus = prep_rls(
1.3400 - Rls {id = "discard_minus", preconds = [], rew_ord = ("dummy_ord",dummy_ord),
1.3401 - erls = e_rls, srls = Erls, calc = [], (*asm_thm = [],*)
1.3402 - rules = [Thm ("real_diff_minus", num_str real_diff_minus),
1.3403 - (*"a - b = a + -1 * b"*)
1.3404 - Thm ("sym_real_mult_minus1",num_str (real_mult_minus1 RS sym))
1.3405 - (*- ?z = "-1 * ?z"*)
1.3406 - ],
1.3407 - scr = Script ((term_of o the o (parse thy))
1.3408 - "empty_script")
1.3409 - }):rls;
1.3410 -(*erls for calculate_Rational; make local with FIXX@ME result:term *term list*)
1.3411 -val powers_erls = prep_rls(
1.3412 - Rls {id = "powers_erls", preconds = [], rew_ord = ("dummy_ord",dummy_ord),
1.3413 - erls = e_rls, srls = Erls, calc = [], (*asm_thm = [],*)
1.3414 - rules = [Calc ("Atools.is'_atom",eval_is_atom "#is_atom_"),
1.3415 - Calc ("Atools.is'_even",eval_is_even "#is_even_"),
1.3416 - Calc ("op <",eval_equ "#less_"),
1.3417 - Thm ("not_false", not_false),
1.3418 - Thm ("not_true", not_true),
1.3419 - Calc ("op +",eval_binop "#add_")
1.3420 - ],
1.3421 - scr = Script ((term_of o the o (parse thy))
1.3422 - "empty_script")
1.3423 - }:rls);
1.3424 -(*.all powers over + distributed; atoms over * collected, other distributed
1.3425 - contains absolute minimum of thms for context in norm_Rational .*)
1.3426 -val powers = prep_rls(
1.3427 - Rls {id = "powers", preconds = [], rew_ord = ("dummy_ord",dummy_ord),
1.3428 - erls = powers_erls, srls = Erls, calc = [], (*asm_thm = [],*)
1.3429 - rules = [Thm ("realpow_multI", num_str realpow_multI),
1.3430 - (*"(r * s) ^^^ n = r ^^^ n * s ^^^ n"*)
1.3431 - Thm ("realpow_pow",num_str realpow_pow),
1.3432 - (*"(a ^^^ b) ^^^ c = a ^^^ (b * c)"*)
1.3433 - Thm ("realpow_oneI",num_str realpow_oneI),
1.3434 - (*"r ^^^ 1 = r"*)
1.3435 - Thm ("realpow_minus_even",num_str realpow_minus_even),
1.3436 - (*"n is_even ==> (- r) ^^^ n = r ^^^ n" ?-->discard_minus?*)
1.3437 - Thm ("realpow_minus_odd",num_str realpow_minus_odd),
1.3438 - (*"Not (n is_even) ==> (- r) ^^^ n = -1 * r ^^^ n"*)
1.3439 -
1.3440 - (*----- collect atoms over * -----*)
1.3441 - Thm ("realpow_two_atom",num_str realpow_two_atom),
1.3442 - (*"r is_atom ==> r * r = r ^^^ 2"*)
1.3443 - Thm ("realpow_plus_1",num_str realpow_plus_1),
1.3444 - (*"r is_atom ==> r * r ^^^ n = r ^^^ (n + 1)"*)
1.3445 - Thm ("realpow_addI_atom",num_str realpow_addI_atom),
1.3446 - (*"r is_atom ==> r ^^^ n * r ^^^ m = r ^^^ (n + m)"*)
1.3447 -
1.3448 - (*----- distribute none-atoms -----*)
1.3449 - Thm ("realpow_def_atom",num_str realpow_def_atom),
1.3450 - (*"[| 1 < n; not(r is_atom) |]==>r ^^^ n = r * r ^^^ (n + -1)"*)
1.3451 - Thm ("realpow_eq_oneI",num_str realpow_eq_oneI),
1.3452 - (*"1 ^^^ n = 1"*)
1.3453 - Calc ("op +",eval_binop "#add_")
1.3454 - ],
1.3455 - scr = Script ((term_of o the o (parse thy))
1.3456 - "empty_script")
1.3457 - }:rls);
1.3458 -(*.contains absolute minimum of thms for context in norm_Rational.*)
1.3459 -val rat_mult_divide = prep_rls(
1.3460 - Rls {id = "rat_mult_divide", preconds = [],
1.3461 - rew_ord = ("dummy_ord",dummy_ord),
1.3462 - erls = e_rls, srls = Erls, calc = [], (*asm_thm = [],*)
1.3463 - rules = [Thm ("rat_mult",num_str rat_mult),
1.3464 - (*(1)"?a / ?b * (?c / ?d) = ?a * ?c / (?b * ?d)"*)
1.3465 - Thm ("real_times_divide1_eq",num_str real_times_divide1_eq),
1.3466 - (*(2)"?a * (?c / ?d) = ?a * ?c / ?d" must be [2],
1.3467 - otherwise inv.to a / b / c = ...*)
1.3468 - Thm ("real_times_divide2_eq",num_str real_times_divide2_eq),
1.3469 - (*"?a / ?b * ?c = ?a * ?c / ?b" order weights x^^^n too much
1.3470 - and does not commute a / b * c ^^^ 2 !*)
1.3471 -
1.3472 - Thm ("real_divide_divide1_eq", real_divide_divide1_eq),
1.3473 - (*"?x / (?y / ?z) = ?x * ?z / ?y"*)
1.3474 - Thm ("real_divide_divide2_eq", real_divide_divide2_eq),
1.3475 - (*"?x / ?y / ?z = ?x / (?y * ?z)"*)
1.3476 - Calc ("HOL.divide" ,eval_cancel "#divide_")
1.3477 - ],
1.3478 - scr = Script ((term_of o the o (parse thy)) "empty_script")
1.3479 - }:rls);
1.3480 -(*.contains absolute minimum of thms for context in norm_Rational.*)
1.3481 -val reduce_0_1_2 = prep_rls(
1.3482 - Rls{id = "reduce_0_1_2", preconds = [], rew_ord = ("dummy_ord", dummy_ord),
1.3483 - erls = e_rls,srls = Erls,calc = [],(*asm_thm = [],*)
1.3484 - rules = [(*Thm ("real_divide_1",num_str real_divide_1),
1.3485 - "?x / 1 = ?x" unnecess.for normalform*)
1.3486 - Thm ("real_mult_1",num_str real_mult_1),
1.3487 - (*"1 * z = z"*)
1.3488 - (*Thm ("real_mult_minus1",num_str real_mult_minus1),
1.3489 - "-1 * z = - z"*)
1.3490 - (*Thm ("real_minus_mult_cancel",num_str real_minus_mult_cancel),
1.3491 - "- ?x * - ?y = ?x * ?y"*)
1.3492 -
1.3493 - Thm ("real_mult_0",num_str real_mult_0),
1.3494 - (*"0 * z = 0"*)
1.3495 - Thm ("real_add_zero_left",num_str real_add_zero_left),
1.3496 - (*"0 + z = z"*)
1.3497 - (*Thm ("real_add_minus",num_str real_add_minus),
1.3498 - "?z + - ?z = 0"*)
1.3499 -
1.3500 - Thm ("sym_real_mult_2",num_str (real_mult_2 RS sym)),
1.3501 - (*"z1 + z1 = 2 * z1"*)
1.3502 - Thm ("real_mult_2_assoc",num_str real_mult_2_assoc),
1.3503 - (*"z1 + (z1 + k) = 2 * z1 + k"*)
1.3504 -
1.3505 - Thm ("real_0_divide",num_str real_0_divide)
1.3506 - (*"0 / ?x = 0"*)
1.3507 - ], scr = EmptyScr}:rls);
1.3508 -
1.3509 -(*erls for calculate_Rational;
1.3510 - make local with FIXX@ME result:term *term list WN0609???SKMG*)
1.3511 -val norm_rat_erls = prep_rls(
1.3512 - Rls {id = "norm_rat_erls", preconds = [], rew_ord = ("dummy_ord",dummy_ord),
1.3513 - erls = e_rls, srls = Erls, calc = [], (*asm_thm = [],*)
1.3514 - rules = [Calc ("Atools.is'_const",eval_const "#is_const_")
1.3515 - ],
1.3516 - scr = Script ((term_of o the o (parse thy))
1.3517 - "empty_script")
1.3518 - }:rls);
1.3519 -(*.consists of rls containing the absolute minimum of thms.*)
1.3520 -(*040209: this version has been used by RL for his equations,
1.3521 -which is now replaced by MGs version below
1.3522 -vvv OLD VERSION !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!*)
1.3523 -val norm_Rational = prep_rls(
1.3524 - Rls {id = "norm_Rational", preconds = [], rew_ord = ("dummy_ord",dummy_ord),
1.3525 - erls = norm_rat_erls, srls = Erls, calc = [], (*asm_thm = [],*)
1.3526 - rules = [(*sequence given by operator precedence*)
1.3527 - Rls_ discard_minus,
1.3528 - Rls_ powers,
1.3529 - Rls_ rat_mult_divide,
1.3530 - Rls_ expand,
1.3531 - Rls_ reduce_0_1_2,
1.3532 - (*^^^^^^^^^ from RL -- not the latest one vvvvvvvvv*)
1.3533 - Rls_ order_add_mult,
1.3534 - Rls_ collect_numerals,
1.3535 - Rls_ add_fractions_p,
1.3536 - Rls_ cancel_p
1.3537 - ],
1.3538 - scr = Script ((term_of o the o (parse thy))
1.3539 - "empty_script")
1.3540 - }:rls);
1.3541 -val norm_Rational_parenthesized = prep_rls(
1.3542 - Seq {id = "norm_Rational_parenthesized", preconds = []:term list,
1.3543 - rew_ord = ("dummy_ord", dummy_ord),
1.3544 - erls = Atools_erls, srls = Erls,
1.3545 - calc = [], (*asm_thm = [],*)
1.3546 - rules = [Rls_ norm_Rational, (*from RL -- not the latest one*)
1.3547 - Rls_ discard_parentheses
1.3548 - ],
1.3549 - scr = EmptyScr
1.3550 - }:rls);
1.3551 -
1.3552 -
1.3553 -(*-------------------18.3.03 --> struct <-----------^^^--*)
1.3554 -
1.3555 -
1.3556 -
1.3557 -theory' := overwritel (!theory', [("Rational.thy",Rational.thy)]);
1.3558 -
1.3559 -
1.3560 -(*WN030318???SK: simplifies all but cancel and common_nominator*)
1.3561 -val simplify_rational =
1.3562 - merge_rls "simplify_rational" expand_binoms
1.3563 - (append_rls "divide" calculate_Rational
1.3564 - [Thm ("real_divide_1",num_str real_divide_1),
1.3565 - (*"?x / 1 = ?x"*)
1.3566 - Thm ("rat_mult",num_str rat_mult),
1.3567 - (*(1)"?a / ?b * (?c / ?d) = ?a * ?c / (?b * ?d)"*)
1.3568 - Thm ("real_times_divide1_eq",num_str real_times_divide1_eq),
1.3569 - (*(2)"?a * (?c / ?d) = ?a * ?c / ?d" must be [2],
1.3570 - otherwise inv.to a / b / c = ...*)
1.3571 - Thm ("real_times_divide2_eq",num_str real_times_divide2_eq),
1.3572 - (*"?a / ?b * ?c = ?a * ?c / ?b"*)
1.3573 - Thm ("add_minus",num_str add_minus),
1.3574 - (*"?a + ?b - ?b = ?a"*)
1.3575 - Thm ("add_minus1",num_str add_minus1),
1.3576 - (*"?a - ?b + ?b = ?a"*)
1.3577 - Thm ("real_divide_minus1",num_str real_divide_minus1)
1.3578 - (*"?x / -1 = - ?x"*)
1.3579 -(*
1.3580 -,
1.3581 - Thm ("",num_str )
1.3582 -*)
1.3583 - ]);
1.3584 -
1.3585 -(*---------vvv-------------MG ab 1.07.2003--------------vvv-----------*)
1.3586 -
1.3587 -(* ------------------------------------------------------------------ *)
1.3588 -(* Simplifier für beliebige Buchterme *)
1.3589 -(* ------------------------------------------------------------------ *)
1.3590 -(*----------------------- norm_Rational_mg ---------------------------*)
1.3591 -(*. description of the simplifier see MG-DA.p.56ff .*)
1.3592 -(* ------------------------------------------------------------------- *)
1.3593 -val common_nominator_p_rls = prep_rls(
1.3594 - Rls {id = "common_nominator_p_rls", preconds = [],
1.3595 - rew_ord = ("dummy_ord",dummy_ord),
1.3596 - erls = e_rls, srls = Erls, calc = [],
1.3597 - rules =
1.3598 - [Rls_ common_nominator_p
1.3599 - (*FIXME.WN0401 ? redesign Rrls - use exhaustively on a term ?
1.3600 - FIXME.WN0510 unnecessary nesting: introduce RRls_ : rls -> rule*)
1.3601 - ],
1.3602 - scr = EmptyScr});
1.3603 -(* ------------------------------------------------------------------- *)
1.3604 -val cancel_p_rls = prep_rls(
1.3605 - Rls {id = "cancel_p_rls", preconds = [],
1.3606 - rew_ord = ("dummy_ord",dummy_ord),
1.3607 - erls = e_rls, srls = Erls, calc = [],
1.3608 - rules =
1.3609 - [Rls_ cancel_p
1.3610 - (*FIXME.WN.0401 ? redesign Rrls - use exhaustively on a term ?*)
1.3611 - ],
1.3612 - scr = EmptyScr});
1.3613 -(* -------------------------------------------------------------------- *)
1.3614 -(*. makes 'normal' fractions; 'is_polyexp' inhibits double fractions;
1.3615 - used in initial part norm_Rational_mg, see example DA-M02-main.p.60.*)
1.3616 -val rat_mult_poly = prep_rls(
1.3617 - Rls {id = "rat_mult_poly", preconds = [],
1.3618 - rew_ord = ("dummy_ord",dummy_ord),
1.3619 - erls = append_rls "e_rls-is_polyexp" e_rls
1.3620 - [Calc ("Poly.is'_polyexp", eval_is_polyexp "")],
1.3621 - srls = Erls, calc = [],
1.3622 - rules =
1.3623 - [Thm ("rat_mult_poly_l",num_str rat_mult_poly_l),
1.3624 - (*"?c is_polyexp ==> ?c * (?a / ?b) = ?c * ?a / ?b"*)
1.3625 - Thm ("rat_mult_poly_r",num_str rat_mult_poly_r)
1.3626 - (*"?c is_polyexp ==> ?a / ?b * ?c = ?a * ?c / ?b"*)
1.3627 - ],
1.3628 - scr = EmptyScr});
1.3629 -(* ------------------------------------------------------------------ *)
1.3630 -(*. makes 'normal' fractions; 'is_polyexp' inhibits double fractions;
1.3631 - used in looping part norm_Rational_rls, see example DA-M02-main.p.60
1.3632 - .. WHERE THE LATTER DOES ALWAYS WORK, BECAUSE erls = e_rls,
1.3633 - I.E. THE RESPECTIVE ASSUMPTION IS STORED AND Thm APPLIED; WN051028
1.3634 - ... WN0609???MG.*)
1.3635 -val rat_mult_div_pow = prep_rls(
1.3636 - Rls {id = "rat_mult_div_pow", preconds = [],
1.3637 - rew_ord = ("dummy_ord",dummy_ord),
1.3638 - erls = e_rls,
1.3639 - (*FIXME.WN051028 append_rls "e_rls-is_polyexp" e_rls
1.3640 - [Calc ("Poly.is'_polyexp", eval_is_polyexp "")],
1.3641 - with this correction ^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^ we get
1.3642 - error "rational.sml.sml: diff.behav. in norm_Rational_mg 29" etc.
1.3643 - thus we decided to go on with this flaw*)
1.3644 - srls = Erls, calc = [],
1.3645 - rules = [Thm ("rat_mult",num_str rat_mult),
1.3646 - (*"?a / ?b * (?c / ?d) = ?a * ?c / (?b * ?d)"*)
1.3647 - Thm ("rat_mult_poly_l",num_str rat_mult_poly_l),
1.3648 - (*"?c is_polyexp ==> ?c * (?a / ?b) = ?c * ?a / ?b"*)
1.3649 - Thm ("rat_mult_poly_r",num_str rat_mult_poly_r),
1.3650 - (*"?c is_polyexp ==> ?a / ?b * ?c = ?a * ?c / ?b"*)
1.3651 -
1.3652 - Thm ("real_divide_divide1_mg", real_divide_divide1_mg),
1.3653 - (*"y ~= 0 ==> (u / v) / (y / z) = (u * z) / (y * v)"*)
1.3654 - Thm ("real_divide_divide1_eq", real_divide_divide1_eq),
1.3655 - (*"?x / (?y / ?z) = ?x * ?z / ?y"*)
1.3656 - Thm ("real_divide_divide2_eq", real_divide_divide2_eq),
1.3657 - (*"?x / ?y / ?z = ?x / (?y * ?z)"*)
1.3658 - Calc ("HOL.divide" ,eval_cancel "#divide_"),
1.3659 -
1.3660 - Thm ("rat_power", num_str rat_power)
1.3661 - (*"(?a / ?b) ^^^ ?n = ?a ^^^ ?n / ?b ^^^ ?n"*)
1.3662 - ],
1.3663 - scr = Script ((term_of o the o (parse thy)) "empty_script")
1.3664 - }:rls);
1.3665 -(* ------------------------------------------------------------------ *)
1.3666 -val rat_reduce_1 = prep_rls(
1.3667 - Rls {id = "rat_reduce_1", preconds = [],
1.3668 - rew_ord = ("dummy_ord",dummy_ord),
1.3669 - erls = e_rls, srls = Erls, calc = [],
1.3670 - rules = [Thm ("real_divide_1",num_str real_divide_1),
1.3671 - (*"?x / 1 = ?x"*)
1.3672 - Thm ("real_mult_1",num_str real_mult_1)
1.3673 - (*"1 * z = z"*)
1.3674 - ],
1.3675 - scr = Script ((term_of o the o (parse thy)) "empty_script")
1.3676 - }:rls);
1.3677 -(* ------------------------------------------------------------------ *)
1.3678 -(*. looping part of norm_Rational(*_mg*) .*)
1.3679 -val norm_Rational_rls = prep_rls(
1.3680 - Rls {id = "norm_Rational_rls", preconds = [],
1.3681 - rew_ord = ("dummy_ord",dummy_ord),
1.3682 - erls = norm_rat_erls, srls = Erls, calc = [],
1.3683 - rules = [Rls_ common_nominator_p_rls,
1.3684 - Rls_ rat_mult_div_pow,
1.3685 - Rls_ make_rat_poly_with_parentheses,
1.3686 - Rls_ cancel_p_rls,(*FIXME:cancel_p does NOT order sometimes*)
1.3687 - Rls_ rat_reduce_1
1.3688 - ],
1.3689 - scr = Script ((term_of o the o (parse thy)) "empty_script")
1.3690 - }:rls);
1.3691 -(* ------------------------------------------------------------------ *)
1.3692 -(*040109 'norm_Rational'(by RL) replaced by 'norm_Rational_mg'(MG)
1.3693 - just be renaming:*)
1.3694 -val norm_Rational(*_mg*) = prep_rls(
1.3695 - Seq {id = "norm_Rational"(*_mg*), preconds = [],
1.3696 - rew_ord = ("dummy_ord",dummy_ord),
1.3697 - erls = norm_rat_erls, srls = Erls, calc = [],
1.3698 - rules = [Rls_ discard_minus_,
1.3699 - Rls_ rat_mult_poly,(* removes double fractions like a/b/c *)
1.3700 - Rls_ make_rat_poly_with_parentheses, (*WN0510 also in(#)below*)
1.3701 - Rls_ cancel_p_rls, (*FIXME.MG:cancel_p does NOT order sometim*)
1.3702 - Rls_ norm_Rational_rls, (* the main rls, looping (#) *)
1.3703 - Rls_ discard_parentheses_ (* mult only *)
1.3704 - ],
1.3705 - scr = Script ((term_of o the o (parse thy)) "empty_script")
1.3706 - }:rls);
1.3707 -(* ------------------------------------------------------------------ *)
1.3708 -
1.3709 -
1.3710 -ruleset' := overwritelthy thy (!ruleset',
1.3711 - [("calculate_Rational", calculate_Rational),
1.3712 - ("calc_rat_erls",calc_rat_erls),
1.3713 - ("rational_erls", rational_erls),
1.3714 - ("cancel_p", cancel_p),
1.3715 - ("cancel", cancel),
1.3716 - ("common_nominator_p", common_nominator_p),
1.3717 - ("common_nominator_p_rls", common_nominator_p_rls),
1.3718 - ("common_nominator" , common_nominator),
1.3719 - ("discard_minus", discard_minus),
1.3720 - ("powers_erls", powers_erls),
1.3721 - ("powers", powers),
1.3722 - ("rat_mult_divide", rat_mult_divide),
1.3723 - ("reduce_0_1_2", reduce_0_1_2),
1.3724 - ("rat_reduce_1", rat_reduce_1),
1.3725 - ("norm_rat_erls", norm_rat_erls),
1.3726 - ("norm_Rational", norm_Rational),
1.3727 - ("norm_Rational_rls", norm_Rational_rls),
1.3728 - ("norm_Rational_parenthesized", norm_Rational_parenthesized),
1.3729 - ("rat_mult_poly", rat_mult_poly),
1.3730 - ("rat_mult_div_pow", rat_mult_div_pow),
1.3731 - ("cancel_p_rls", cancel_p_rls)
1.3732 - ]);
1.3733 -
1.3734 -calclist':= overwritel (!calclist',
1.3735 - [("is_expanded", ("Rational.is'_expanded", eval_is_expanded ""))
1.3736 - ]);
1.3737 -
1.3738 -(** problems **)
1.3739 -
1.3740 -store_pbt
1.3741 - (prep_pbt Rational.thy "pbl_simp_rat" [] e_pblID
1.3742 - (["rational","simplification"],
1.3743 - [("#Given" ,["term t_"]),
1.3744 - ("#Where" ,["t_ is_ratpolyexp"]),
1.3745 - ("#Find" ,["normalform n_"])
1.3746 - ],
1.3747 - append_rls "e_rls" e_rls [(*for preds in where_*)],
1.3748 - SOME "Simplify t_",
1.3749 - [["simplification","of_rationals"]]));
1.3750 -
1.3751 -(** methods **)
1.3752 -
1.3753 -(*WN061025 this methods script is copied from (auto-generated) script
1.3754 - of norm_Rational in order to ease repair on inform*)
1.3755 -store_met
1.3756 - (prep_met Rational.thy "met_simp_rat" [] e_metID
1.3757 - (["simplification","of_rationals"],
1.3758 - [("#Given" ,["term t_"]),
1.3759 - ("#Where" ,["t_ is_ratpolyexp"]),
1.3760 - ("#Find" ,["normalform n_"])
1.3761 - ],
1.3762 - {rew_ord'="tless_true",
1.3763 - rls' = e_rls,
1.3764 - calc = [], srls = e_rls,
1.3765 - prls = append_rls "simplification_of_rationals_prls" e_rls
1.3766 - [(*for preds in where_*)
1.3767 - Calc ("Rational.is'_ratpolyexp",
1.3768 - eval_is_ratpolyexp "")],
1.3769 - crls = e_rls, nrls = norm_Rational_rls},
1.3770 -"Script SimplifyScript (t_::real) = \
1.3771 -\ ((Try (Rewrite_Set discard_minus_ False) @@ \
1.3772 -\ Try (Rewrite_Set rat_mult_poly False) @@ \
1.3773 -\ Try (Rewrite_Set make_rat_poly_with_parentheses False) @@ \
1.3774 -\ Try (Rewrite_Set cancel_p_rls False) @@ \
1.3775 -\ (Repeat \
1.3776 -\ ((Try (Rewrite_Set common_nominator_p_rls False) @@ \
1.3777 -\ Try (Rewrite_Set rat_mult_div_pow False) @@ \
1.3778 -\ Try (Rewrite_Set make_rat_poly_with_parentheses False) @@\
1.3779 -\ Try (Rewrite_Set cancel_p_rls False) @@ \
1.3780 -\ Try (Rewrite_Set rat_reduce_1 False)))) @@ \
1.3781 -\ Try (Rewrite_Set discard_parentheses_ False)) \
1.3782 -\ t_)"
1.3783 - ));
1.3784 -
1.3785 -(* use"../IsacKnowledge/Rational.ML";
1.3786 - use"IsacKnowledge/Rational.ML";
1.3787 - use"Rational.ML";
1.3788 - *)
1.3789 -