(* tools for integration over the reals

   author: Walther Neuper 050905, 08:51

   (c) due to copyright terms

use"IsacKnowledge/Integrate.ML";

use"Integrate.ML";

*)

(** interface isabelle -- isac **)

theory' := overwritel (!theory', [("Integrate.thy",Integrate.thy)]);

theorem' := overwritel (!theorem',

[("integral_const",num_str integral_const),

 ("integral_var",num_str integral_var),

 ("integral_add",num_str integral_add),

 ("integral_mult",num_str integral_mult),

 ("call_for_new_c",num_str call_for_new_c),

 ("integral_pow",num_str integral_pow)

]);

(** eval functions **)

val c = Free ("c", HOLogic.realT);

(*.create a new unique variable 'c..' in a term; for use by Calc in a rls;

   an alternative to do this would be '(Try (Calculate new_c_) (new_c es__))'

   in the script; this will be possible if currying doesnt take the value

   from a variable, but the value '(new_c es__)' itself.*)

fun new_c term = 

    let fun selc var = 

	    case (explode o id_of) var of

		"c"::[] => true

	      |	"c"::"_"::is => (case (int_of_str o implode) is of

				     Some _ => true

				   | None => false)

              | _ => false;

	fun get_coeff c = case (explode o id_of) c of

	      		      "c"::"_"::is => (the o int_of_str o implode) is

			    | _ => 0;

        val cs = filter selc (vars term);

    in 

	case cs of

	    [] => c

	  | [c] => Free ("c_2", HOLogic.realT)

	  | cs => 

	    let val max_coeff = maxl (map get_coeff cs)

	    in Free ("c_"^string_of_int (max_coeff + 1), HOLogic.realT) end

    end;

(*("new_c", ("Integrate.new'_c", eval_new_c "#new_c_"))*)

fun eval_new_c _ _ (p as (Const ("Integrate.new'_c",_) $ t)) _ =

     Some ((term2str p) ^ " = " ^ term2str (new_c p),

	  Trueprop $ (mk_equality (p, new_c p)))

  | eval_new_c _ _ _ _ = None;

(*("is_f_x", ("Integrate.is'_f'_x", eval_is_f_x "is_f_x_"))*)

fun eval_is_f_x _ _(p as (Const ("Integrate.is'_f'_x", _)

					   $ arg)) _ =

    if is_f_x arg

    then Some ((term2str p) ^ " = True",

	       Trueprop $ (mk_equality (p, HOLogic.true_const)))

    else Some ((term2str p) ^ " = False",

	       Trueprop $ (mk_equality (p, HOLogic.false_const)))

  | eval_is_f_x _ _ _ _ = None;

calclist':= overwritel (!calclist', 

   [("new_c", ("Integrate.new'_c", eval_new_c "new_c_")),

    ("is_f_x", ("Integrate.is'_f'_x", eval_is_f_x "is_f_idextifier_"))

    ]);

(** rewrite order 'ord_simplify_Integral' **)

(* order wrt. several linear (i.e. without exponents) variables "c","c_2",..

   which leaves the monomials containing c, c_2,... at the end of an Integral

   and puts the c, c_2,... rightmost within a monomial.

   WN050906 this is a quick and dirty adaption of ord_make_polynomial_in,

   which was most adequate, because it uses size_of_term*)

(**)

local (*. for simplify_Integral .*)

(**)

open Term;  (* for type order = EQUAL | LESS | GREATER *)

fun pr_ord EQUAL = "EQUAL"

  | pr_ord LESS  = "LESS"

  | pr_ord GREATER = "GREATER";

fun dest_hd' (Const (a, T)) = (((a, 0), T), 0)

  | dest_hd' (Free (ccc, T)) =

    (case explode ccc of

	"c"::[] => ((("|||||||||||||||||||||", 0), T), 1)(*greatest string WN*)

      | "c"::"_"::_ => ((("|||||||||||||||||||||", 0), T), 1)

      | _ => (((ccc, 0), T), 1))

  | dest_hd' (Var v) = (v, 2)

  | dest_hd' (Bound i) = ((("", i), dummyT), 3)

  | dest_hd' (Abs (_, T, _)) = ((("", 0), T), 4);

fun size_of_term' (Free (ccc, _)) =

    (case explode ccc of

	"c"::[] => 1000

      | "c"::"_"::is => 1000 * ((str2int o implode) is)

      | _ => 1)

  | size_of_term' (Abs (_,_,body)) = 1 + size_of_term' body

  | size_of_term' (f$t) = size_of_term' f  +  size_of_term' t

  | size_of_term' _ = 1;

fun term_ord' pr thy (Abs (_, T, t), Abs(_, U, u)) =       (* ~ term.ML *)

      (case term_ord' pr thy (t, u) of EQUAL => typ_ord (T, U) | ord => ord)

  | term_ord' pr thy (t, u) =

      (if pr then 

	 let

	   val (f, ts) = strip_comb t and (g, us) = strip_comb u;

	   val _=writeln("t= f@ts= \""^

	      ((string_of_cterm o cterm_of (sign_of thy)) f)^"\" @ \"["^

	      (commas(map(string_of_cterm o cterm_of(sign_of thy)) ts))^"]\"");

	   val _=writeln("u= g@us= \""^

	      ((string_of_cterm o cterm_of (sign_of thy)) g)^"\" @ \"["^

	      (commas(map(string_of_cterm o cterm_of(sign_of thy)) us))^"]\"");

	   val _=writeln("size_of_term(t,u)= ("^

	      (string_of_int(size_of_term' t))^", "^

	      (string_of_int(size_of_term' u))^")");

	   val _=writeln("hd_ord(f,g)      = "^((pr_ord o hd_ord)(f,g)));

	   val _=writeln("terms_ord(ts,us) = "^

			   ((pr_ord o terms_ord str false)(ts,us)));

	   val _=writeln("-------");

	 in () end

       else ();

	 case int_ord (size_of_term' t, size_of_term' u) of

	   EQUAL =>

	     let val (f, ts) = strip_comb t and (g, us) = strip_comb u in

	       (case hd_ord (f, g) of EQUAL => (terms_ord str pr) (ts, us) 

	     | ord => ord)

	     end

	 | ord => ord)

and hd_ord (f, g) =                                        (* ~ term.ML *)

  prod_ord (prod_ord indexname_ord typ_ord) int_ord (dest_hd' f, 

						     dest_hd' g)

and terms_ord str pr (ts, us) = 

    list_ord (term_ord' pr (assoc_thy "Isac.thy"))(ts, us);

(**)

in

(**)

(*WN0510 for preliminary use in eval_order_system, see case-study mat-eng.tex

fun ord_simplify_Integral_rev (pr:bool) thy subst tu = 

    (term_ord' pr thy (Library.swap tu) = LESS);*)

(*for the rls's*)

fun ord_simplify_Integral (pr:bool) thy subst tu = 

    (term_ord' pr thy tu = LESS);

(**)

end;

(**)

rew_ord' := overwritel (!rew_ord',

[("ord_simplify_Integral", ord_simplify_Integral false thy)

 ]);

(** rulesets **)

(*.rulesets for integration.*)

val integration_rules = 

    prep_rls (Rls {id="integration_rules", preconds = [], 

		   rew_ord = ("termlessI",termlessI), 

		   erls = Rls {id="conditions_in_integration_rules", 

			       preconds = [], 

			       rew_ord = ("termlessI",termlessI), 

			       erls = Erls, 

			       srls = Erls, calc = [],

			       rules = [(*for rewriting conditions in Thm's*)

					Calc ("Atools.occurs'_in", 

					      eval_occurs_in "#occurs_in_"),

					Thm ("not_true",num_str not_true),

					Thm ("not_false",not_false)

					],

			       scr = EmptyScr}, 

		   srls = Erls, calc = [],

		   rules = [

			    Thm ("integral_const",num_str integral_const),

			    Thm ("integral_var",num_str integral_var),

			    Thm ("integral_add",num_str integral_add),

			    Thm ("integral_mult",num_str integral_mult),

			    Thm ("integral_pow",num_str integral_pow),

			    Calc ("op +", eval_binop "#add_")(*for n+1*)

			    ],

		   scr = EmptyScr});

val add_new_c = 

    prep_rls (Seq {id="add_new_c", preconds = [], 

		   rew_ord = ("termlessI",termlessI), 

		   erls = Rls {id="conditions_in_add_new_c", 

			       preconds = [], 

			       rew_ord = ("termlessI",termlessI), 

			       erls = Erls, 

			       srls = Erls, calc = [],

			       rules = [Calc ("Tools.matches", eval_matches""),

					Calc ("Integrate.is'_f'_x", 

					 eval_is_f_x "is_f_x_"),

					Thm ("not_true",num_str not_true),

					Thm ("not_false",num_str not_false)

					],

			       scr = EmptyScr}, 

		   srls = Erls, calc = [],

		   rules = [ Thm ("call_for_new_c", num_str call_for_new_c),

			     Calc("Integrate.new'_c", eval_new_c "new_c_")

			    ],

		   scr = EmptyScr});

(*.rulesets for simplifying Integrals.*)

(*'' is adapted from '' by just replacing the rew_ord*)

val order_add_mult_Integral = 

  Rls{id = "order_add_mult_Integral", preconds = [], 

      rew_ord = ("ord_simplify_Integral",

		 ord_simplify_Integral false Integrate.thy),

      erls = e_rls,srls = Erls, calc = [],

      rules = [Thm ("real_mult_commute",num_str real_mult_commute),

	       (* z * w = w * z *)

	       Thm ("real_mult_left_commute",num_str real_mult_left_commute),

	       (*z1.0 * (z2.0 * z3.0) = z2.0 * (z1.0 * z3.0)*)

	       Thm ("real_mult_assoc",num_str real_mult_assoc),		

	       (*z1.0 * z2.0 * z3.0 = z1.0 * (z2.0 * z3.0)*)

	       Thm ("real_add_commute",num_str real_add_commute),	

	       (*z + w = w + z*)

	       Thm ("real_add_left_commute",num_str real_add_left_commute),

	       (*x + (y + z) = y + (x + z)*)

	       Thm ("real_add_assoc",num_str real_add_assoc)	               

	       (*z1.0 + z2.0 + z3.0 = z1.0 + (z2.0 + z3.0)*)

	       ], 

      scr = EmptyScr}:rls;

(*'norm_Rational_no_add_fractions' is adapted from 'norm_Rational' by

  #1 using 'ord_simplify_Integral' in 'order_add_mult_Integral'

  #2 NOT using common_nominator_p

  *)

val norm_Rational_no_add_fractions = prep_rls(

  Rls {id = "norm_Rational_no_add_fractions", preconds = [], 

       rew_ord = ("dummy_ord",dummy_ord), 

       erls = norm_rat_erls, srls = Erls, calc = [],

       rules = [(*sequence given by operator precedence*)

		Rls_ discard_minus,

		Rls_ powers,

		Rls_ rat_mult_divide,

		Rls_ expand,

		Rls_ reduce_0_1_2,

		Rls_ (*order_add_mult #1*) order_add_mult_Integral,

		Rls_ collect_numerals,

		(*Rls_ add_fractions_p, #2*)

		Rls_ cancel_p

		],

       scr = Script ((term_of o the o (parse thy)) 

			 "empty_script")

       }:rls);

(*'simplify_Integral' is adapted from 'norm_Rational_parenthesized'

  this simplifier orders 'c'..'c_123' right most and 

  does NOT 'common_nominator' for fractions*)

val simplify_Integral = prep_rls(

  Seq {id = "simplify_Integral", preconds = []:term list, 

       rew_ord = ("dummy_ord", dummy_ord),

      erls = Atools_erls, srls = Erls, calc = [],

      rules = [Rls_ norm_Rational_no_add_fractions,

	       Rls_ discard_parentheses

	       ],

      scr = EmptyScr

      }:rls);      

(*rls (ab-)used preliminarily for 'EqSystem's in 'Biegelininie'*)

val simplify_Integral_parenthesized = prep_rls(

  Seq {id = "simplify_Integral_parenthesized", preconds = []:term list, 

       rew_ord = ("dummy_ord", dummy_ord),

      erls = Atools_erls, srls = Erls, calc = [],

      rules = [Rls_ norm_Rational_no_add_fractions(*,

	       Rls_ discard_parentheses*)

	       ],

      scr = EmptyScr

      }:rls);      

val integration = 

    prep_rls (Seq {id="integration", preconds = [], 

		   rew_ord = ("termlessI",termlessI), 

		   erls = Rls {id="conditions_in_integration", 

			       preconds = [], 

			       rew_ord = ("termlessI",termlessI), 

			       erls = Erls, 

			       srls = Erls, calc = [],

			       rules = [],

			       scr = EmptyScr}, 

		   srls = Erls, calc = [],

		   rules = [ Rls_ integration_rules,

			     Rls_ add_new_c,

			     Rls_ simplify_Integral

			    ],

		   scr = EmptyScr});

ruleset' := 

overwritel (!ruleset', 

	    [("integration_rules", integration_rules),

	     ("add_new_c", add_new_c),

	     ("order_add_mult_Integral", order_add_mult_Integral),

	     ("norm_Rational_no_add_fractions", 

	      norm_Rational_no_add_fractions),

	     ("simplify_Integral", simplify_Integral),

	     ("integration", integration)]);

(** problems **)

store_pbt

 (prep_pbt Integrate.thy

 (["integrate","function"],

  [("#Given" ,["functionTerm f_", "integrateBy v_"]),

   ("#Find"  ,["antiDerivative F_"])

  ],

  append_rls "e_rls" e_rls [(*for preds in where_*)], 

  Some "Integrate (f_, v_)", 

  [["Diff","integration"]]));

store_pbt

 (prep_pbt Integrate.thy

 (["named","integrate","function"],

  [("#Given" ,["functionTerm f_", "integrateBy v_"]),

   ("#Find"  ,["antiDerivativeName F_"])

  ],

  append_rls "e_rls" e_rls [(*for preds in where_*)], 

  Some "Integrate (f_, v_)", 

  [["Diff","integration","named"]]));

(** methods **)

store_met

    (prep_met Integrate.thy

	      (["Diff","integration"],

	       [("#Given" ,["functionTerm f_", "integrateBy v_"]),

		("#Find"  ,["antiDerivative F_"])

		],

	       {rew_ord'="tless_true", rls'=Atools_erls, calc = [], 

		srls = e_rls, 

		prls=e_rls,

	     crls = Atools_erls, nrls = e_rls},

"Script IntegrationScript (f_::real) (v_::real) =                \

\  (let t_ = Take (Integral f_ D v_)                             \

\   in (Rewrite_Set_Inst [(bdv,v_)] integration False) (t_::real))"

));

store_met

    (prep_met Integrate.thy

	      (["Diff","integration","named"],

	       [("#Given" ,["functionTerm f_", "integrateBy v_"]),

		("#Find"  ,["antiDerivativeName F_"])

		],

	       {rew_ord'="tless_true", rls'=Atools_erls, calc = [], 

		srls = e_rls, 

		prls=e_rls,

		crls = Atools_erls, nrls = e_rls},

"Script NamedIntegrationScript (f_::real) (v_::real) (F_::real=>real) = \

\  (let t_ = Take (F_ v_ = Integral f_ D v_)                         \

\   in (Rewrite_Set_Inst [(bdv,v_)] integration False) t_)"

 ));