src/Tools/isac/Knowledge/Integrate.thy
author Walther Neuper <neuper@ist.tugraz.at>
Tue, 31 Aug 2010 16:00:13 +0200
branchisac-update-Isa09-2
changeset 37966 78938fc8e022
parent 37965 9c11005c33b8
child 37967 bd4f7a35e892
permissions -rw-r--r--
num_str --> num_str @{thm
     1 (* integration over the reals
     2    author: Walther Neuper
     3    050814, 08:51
     4    (c) due to copyright terms
     5 *)
     6 
     7 theory Integrate imports Diff begin
     8 
     9 consts
    10 
    11   Integral            :: "[real, real]=> real" ("Integral _ D _" 91)
    12 (*new'_c	      :: "real => real"        ("new'_c _" 66)*)
    13   is'_f'_x            :: "real => bool"        ("_ is'_f'_x" 10)
    14 
    15   (*descriptions in the related problems*)
    16   integrateBy         :: real => una
    17   antiDerivative      :: real => una
    18   antiDerivativeName  :: (real => real) => una
    19 
    20   (*the CAS-command, eg. "Integrate (2*x^^^3, x)"*)
    21   Integrate           :: "[real * real] => real"
    22 
    23   (*Script-names*)
    24   IntegrationScript      :: "[real,real,  real] => real"
    25                   ("((Script IntegrationScript (_ _ =))// (_))" 9)
    26   NamedIntegrationScript :: "[real,real, real=>real,  bool] => bool"
    27                   ("((Script NamedIntegrationScript (_ _ _=))// (_))" 9)
    28 
    29 axioms 
    30 (*stated as axioms, todo: prove as theorems
    31   'bdv' is a constant handled on the meta-level 
    32    specifically as a 'bound variable'            *)
    33 
    34   integral_const    "Not (bdv occurs_in u) ==> Integral u D bdv = u * bdv"
    35   integral_var      "Integral bdv D bdv = bdv ^^^ 2 / 2"
    36 
    37   integral_add      "Integral (u + v) D bdv =  
    38 		     (Integral u D bdv) + (Integral v D bdv)"
    39   integral_mult     "[| Not (bdv occurs_in u); bdv occurs_in v |] ==>  
    40 		     Integral (u * v) D bdv = u * (Integral v D bdv)"
    41 (*WN080222: this goes into sub-terms, too ...
    42   call_for_new_c    "[| Not (matches (u + new_c v) a); Not (a is_f_x) |] ==>  
    43 		     a = a + new_c a"
    44 *)
    45   integral_pow      "Integral bdv ^^^ n D bdv = bdv ^^^ (n+1) / (n + 1)"
    46 
    47 ML {*
    48 (** eval functions **)
    49 
    50 val c = Free ("c", HOLogic.realT);
    51 (*.create a new unique variable 'c..' in a term; for use by Calc in a rls;
    52    an alternative to do this would be '(Try (Calculate new_c_) (new_c es__))'
    53    in the script; this will be possible if currying doesnt take the value
    54    from a variable, but the value '(new_c es__)' itself.*)
    55 fun new_c term = 
    56     let fun selc var = 
    57 	    case (explode o id_of) var of
    58 		"c"::[] => true
    59 	      |	"c"::"_"::is => (case (int_of_str o implode) is of
    60 				     SOME _ => true
    61 				   | NONE => false)
    62               | _ => false;
    63 	fun get_coeff c = case (explode o id_of) c of
    64 	      		      "c"::"_"::is => (the o int_of_str o implode) is
    65 			    | _ => 0;
    66         val cs = filter selc (vars term);
    67     in 
    68 	case cs of
    69 	    [] => c
    70 	  | [c] => Free ("c_2", HOLogic.realT)
    71 	  | cs => 
    72 	    let val max_coeff = maxl (map get_coeff cs)
    73 	    in Free ("c_"^string_of_int (max_coeff + 1), HOLogic.realT) end
    74     end;
    75 
    76 (*WN080222
    77 (*("new_c", ("Integrate.new'_c", eval_new_c "#new_c_"))*)
    78 fun eval_new_c _ _ (p as (Const ("Integrate.new'_c",_) $ t)) _ =
    79      SOME ((term2str p) ^ " = " ^ term2str (new_c p),
    80 	  Trueprop $ (mk_equality (p, new_c p)))
    81   | eval_new_c _ _ _ _ = NONE;
    82 *)
    83 
    84 (*WN080222:*)
    85 (*("add_new_c", ("Integrate.add'_new'_c", eval_add_new_c "#add_new_c_"))
    86   add a new c to a term or a fun-equation;
    87   this is _not in_ the term, because only applied to _whole_ term*)
    88 fun eval_add_new_c (_:string) "Integrate.add'_new'_c" p (_:theory) =
    89     let val p' = case p of
    90 		     Const ("op =", T) $ lh $ rh => 
    91 		     Const ("op =", T) $ lh $ mk_add rh (new_c rh)
    92 		   | p => mk_add p (new_c p)
    93     in SOME ((term2str p) ^ " = " ^ term2str p',
    94 	  Trueprop $ (mk_equality (p, p')))
    95     end
    96   | eval_add_new_c _ _ _ _ = NONE;
    97 
    98 
    99 (*("is_f_x", ("Integrate.is'_f'_x", eval_is_f_x "is_f_x_"))*)
   100 fun eval_is_f_x _ _(p as (Const ("Integrate.is'_f'_x", _)
   101 					   $ arg)) _ =
   102     if is_f_x arg
   103     then SOME ((term2str p) ^ " = True",
   104 	       Trueprop $ (mk_equality (p, HOLogic.true_const)))
   105     else SOME ((term2str p) ^ " = False",
   106 	       Trueprop $ (mk_equality (p, HOLogic.false_const)))
   107   | eval_is_f_x _ _ _ _ = NONE;
   108 
   109 calclist':= overwritel (!calclist', 
   110    [(*("new_c", ("Integrate.new'_c", eval_new_c "new_c_")),*)
   111     ("add_new_c", ("Integrate.add'_new'_c", eval_add_new_c "add_new_c_")),
   112     ("is_f_x", ("Integrate.is'_f'_x", eval_is_f_x "is_f_idextifier_"))
   113     ]);
   114 
   115 
   116 (** rulesets **)
   117 
   118 (*.rulesets for integration.*)
   119 val integration_rules = 
   120     Rls {id="integration_rules", preconds = [], 
   121 	 rew_ord = ("termlessI",termlessI), 
   122 	 erls = Rls {id="conditions_in_integration_rules", 
   123 		     preconds = [], 
   124 		     rew_ord = ("termlessI",termlessI), 
   125 		     erls = Erls, 
   126 		     srls = Erls, calc = [],
   127 		     rules = [(*for rewriting conditions in Thm's*)
   128 			      Calc ("Atools.occurs'_in", 
   129 				    eval_occurs_in "#occurs_in_"),
   130 			      Thm ("not_true",num_str @{not_true),
   131 			      Thm ("not_false",not_false)
   132 			      ],
   133 		     scr = EmptyScr}, 
   134 	 srls = Erls, calc = [],
   135 	 rules = [
   136 		  Thm ("integral_const",num_str @{integral_const),
   137 		  Thm ("integral_var",num_str @{integral_var),
   138 		  Thm ("integral_add",num_str @{integral_add),
   139 		  Thm ("integral_mult",num_str @{integral_mult),
   140 		  Thm ("integral_pow",num_str @{integral_pow),
   141 		  Calc ("op +", eval_binop "#add_")(*for n+1*)
   142 		  ],
   143 	 scr = EmptyScr};
   144 val add_new_c = 
   145     Seq {id="add_new_c", preconds = [], 
   146 	 rew_ord = ("termlessI",termlessI), 
   147 	 erls = Rls {id="conditions_in_add_new_c", 
   148 		     preconds = [], 
   149 		     rew_ord = ("termlessI",termlessI), 
   150 		     erls = Erls, 
   151 		     srls = Erls, calc = [],
   152 		     rules = [Calc ("Tools.matches", eval_matches""),
   153 			      Calc ("Integrate.is'_f'_x", 
   154 				    eval_is_f_x "is_f_x_"),
   155 			      Thm ("not_true",num_str @{not_true),
   156 			      Thm ("not_false",num_str @{not_false)
   157 			      ],
   158 		     scr = EmptyScr}, 
   159 	 srls = Erls, calc = [],
   160 	 rules = [ (*Thm ("call_for_new_c", num_str @{call_for_new_c),*)
   161 		   Cal1 ("Integrate.add'_new'_c", eval_add_new_c "new_c_")
   162 		   ],
   163 	 scr = EmptyScr};
   164 
   165 (*.rulesets for simplifying Integrals.*)
   166 
   167 (*.for simplify_Integral adapted from 'norm_Rational_rls'.*)
   168 val norm_Rational_rls_noadd_fractions = 
   169 Rls {id = "norm_Rational_rls_noadd_fractions", preconds = [], 
   170      rew_ord = ("dummy_ord",dummy_ord), 
   171      erls = norm_rat_erls, srls = Erls, calc = [],
   172      rules = [(*Rls_ common_nominator_p_rls,!!!*)
   173 	      Rls_ (*rat_mult_div_pow original corrected WN051028*)
   174 		  (Rls {id = "rat_mult_div_pow", preconds = [], 
   175 		       rew_ord = ("dummy_ord",dummy_ord), 
   176 		       erls = (*FIXME.WN051028 e_rls,*)
   177 		       append_rls "e_rls-is_polyexp" e_rls
   178 				  [Calc ("Poly.is'_polyexp", 
   179 					 eval_is_polyexp "")],
   180 				  srls = Erls, calc = [],
   181 				  rules = [Thm ("rat_mult",num_str @{rat_mult),
   182 	       (*"?a / ?b * (?c / ?d) = ?a * ?c / (?b * ?d)"*)
   183 	       Thm ("rat_mult_poly_l",num_str @{rat_mult_poly_l),
   184 	       (*"?c is_polyexp ==> ?c * (?a / ?b) = ?c * ?a / ?b"*)
   185 	       Thm ("rat_mult_poly_r",num_str @{rat_mult_poly_r),
   186 	       (*"?c is_polyexp ==> ?a / ?b * ?c = ?a * ?c / ?b"*)
   187 
   188 	       Thm ("real_divide_divide1_mg", real_divide_divide1_mg),
   189 	       (*"y ~= 0 ==> (u / v) / (y / z) = (u * z) / (y * v)"*)
   190 	       Thm ("divide_divide_eq_right", real_divide_divide1_eq),
   191 	       (*"?x / (?y / ?z) = ?x * ?z / ?y"*)
   192 	       Thm ("divide_divide_eq_left", real_divide_divide2_eq),
   193 	       (*"?x / ?y / ?z = ?x / (?y * ?z)"*)
   194 	       Calc ("HOL.divide"  ,eval_cancel "#divide_"),
   195 	      
   196 	       Thm ("rat_power", num_str @{rat_power)
   197 		(*"(?a / ?b) ^^^ ?n = ?a ^^^ ?n / ?b ^^^ ?n"*)
   198 	       ],
   199       scr = Script ((term_of o the o (parse thy)) "empty_script")
   200       }),
   201 		Rls_ make_rat_poly_with_parentheses,
   202 		Rls_ cancel_p_rls,(*FIXME:cancel_p does NOT order sometimes*)
   203 		Rls_ rat_reduce_1
   204 		],
   205        scr = Script ((term_of o the o (parse thy)) "empty_script")
   206        }:rls;
   207 
   208 (*.for simplify_Integral adapted from 'norm_Rational'.*)
   209 val norm_Rational_noadd_fractions = 
   210    Seq {id = "norm_Rational_noadd_fractions", preconds = [], 
   211        rew_ord = ("dummy_ord",dummy_ord), 
   212        erls = norm_rat_erls, srls = Erls, calc = [],
   213        rules = [Rls_ discard_minus_,
   214 		Rls_ rat_mult_poly,(* removes double fractions like a/b/c    *)
   215 		Rls_ make_rat_poly_with_parentheses, (*WN0510 also in(#)below*)
   216 		Rls_ cancel_p_rls, (*FIXME.MG:cancel_p does NOT order sometim*)
   217 		Rls_ norm_Rational_rls_noadd_fractions,(* the main rls (#)   *)
   218 		Rls_ discard_parentheses_ (* mult only                       *)
   219 		],
   220        scr = Script ((term_of o the o (parse thy)) "empty_script")
   221        }:rls;
   222 
   223 (*.simplify terms before and after Integration such that  
   224    ..a.x^2/2 + b.x^3/3.. is made to ..a/2.x^2 + b/3.x^3.. (and NO
   225    common denominator as done by norm_Rational or make_ratpoly_in.
   226    This is a copy from 'make_ratpoly_in' with respective reduction of rules and
   227    *1* expand the term, ie. distribute * and / over +
   228 .*)
   229 val separate_bdv2 =
   230     append_rls "separate_bdv2"
   231 	       collect_bdv
   232 	       [Thm ("separate_bdv", num_str @{separate_bdv),
   233 		(*"?a * ?bdv / ?b = ?a / ?b * ?bdv"*)
   234 		Thm ("separate_bdv_n", num_str @{separate_bdv_n),
   235 		Thm ("separate_1_bdv", num_str @{separate_1_bdv),
   236 		(*"?bdv / ?b = (1 / ?b) * ?bdv"*)
   237 		Thm ("separate_1_bdv_n", num_str @{separate_1_bdv_n)(*,
   238 			  (*"?bdv ^^^ ?n / ?b = 1 / ?b * ?bdv ^^^ ?n"*)
   239 			  *****Thm ("nadd_divide_distrib", 
   240 			  *****num_str @{thm nadd_divide_distrib})
   241 			  (*"(?x + ?y) / ?z = ?x / ?z + ?y / ?z"*)----------*)
   242 		];
   243 val simplify_Integral = 
   244   Seq {id = "simplify_Integral", preconds = []:term list, 
   245        rew_ord = ("dummy_ord", dummy_ord),
   246       erls = Atools_erls, srls = Erls,
   247       calc = [], (*asm_thm = [],*)
   248       rules = [Thm ("left_distrib",num_str @{thm left_distrib}),
   249  	       (*"(?z1.0 + ?z2.0) * ?w = ?z1.0 * ?w + ?z2.0 * ?w"*)
   250 	       Thm ("nadd_divide_distrib",num_str @{thm nadd_divide_distrib}),
   251  	       (*"(?x + ?y) / ?z = ?x / ?z + ?y / ?z"*)
   252 	       (*^^^^^ *1* ^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^*)
   253 	       Rls_ norm_Rational_noadd_fractions,
   254 	       Rls_ order_add_mult_in,
   255 	       Rls_ discard_parentheses,
   256 	       (*Rls_ collect_bdv, from make_polynomial_in*)
   257 	       Rls_ separate_bdv2,
   258 	       Calc ("HOL.divide"  ,eval_cancel "#divide_")
   259 	       ],
   260       scr = EmptyScr}:rls;      
   261 
   262 
   263 (*simplify terms before and after Integration such that  
   264    ..a.x^2/2 + b.x^3/3.. is made to ..a/2.x^2 + b/3.x^3.. (and NO
   265    common denominator as done by norm_Rational or make_ratpoly_in.
   266    This is a copy from 'make_polynomial_in' with insertions from 
   267    'make_ratpoly_in' 
   268 THIS IS KEPT FOR COMPARISON ............................................   
   269 * val simplify_Integral = prep_rls(
   270 *   Seq {id = "", preconds = []:term list, 
   271 *        rew_ord = ("dummy_ord", dummy_ord),
   272 *       erls = Atools_erls, srls = Erls,
   273 *       calc = [], (*asm_thm = [],*)
   274 *       rules = [Rls_ expand_poly,
   275 * 	       Rls_ order_add_mult_in,
   276 * 	       Rls_ simplify_power,
   277 * 	       Rls_ collect_numerals,
   278 * 	       Rls_ reduce_012,
   279 * 	       Thm ("realpow_oneI",num_str @{realpow_oneI),
   280 * 	       Rls_ discard_parentheses,
   281 * 	       Rls_ collect_bdv,
   282 * 	       (*below inserted from 'make_ratpoly_in'*)
   283 * 	       Rls_ (append_rls "separate_bdv"
   284 * 			 collect_bdv
   285 * 			 [Thm ("separate_bdv", num_str @{separate_bdv),
   286 * 			  (*"?a * ?bdv / ?b = ?a / ?b * ?bdv"*)
   287 * 			  Thm ("separate_bdv_n", num_str @{separate_bdv_n),
   288 * 			  Thm ("separate_1_bdv", num_str @{separate_1_bdv),
   289 * 			  (*"?bdv / ?b = (1 / ?b) * ?bdv"*)
   290 * 			  Thm ("separate_1_bdv_n", num_str @{separate_1_bdv_n)(*,
   291 * 			  (*"?bdv ^^^ ?n / ?b = 1 / ?b * ?bdv ^^^ ?n"*)
   292 * 			  Thm ("nadd_divide_distrib", 
   293 * 				 num_str @{thm nadd_divide_distrib})
   294 * 			   (*"(?x + ?y) / ?z = ?x / ?z + ?y / ?z"*)*)
   295 * 			  ]),
   296 * 	       Calc ("HOL.divide"  ,eval_cancel "#divide_")
   297 * 	       ],
   298 *       scr = EmptyScr
   299 *       }:rls); 
   300 .......................................................................*)
   301 
   302 val integration = 
   303     Seq {id="integration", preconds = [], 
   304 	 rew_ord = ("termlessI",termlessI), 
   305 	 erls = Rls {id="conditions_in_integration", 
   306 		     preconds = [], 
   307 		     rew_ord = ("termlessI",termlessI), 
   308 		     erls = Erls, 
   309 		     srls = Erls, calc = [],
   310 		     rules = [],
   311 		     scr = EmptyScr}, 
   312 	 srls = Erls, calc = [],
   313 	 rules = [ Rls_ integration_rules,
   314 		   Rls_ add_new_c,
   315 		   Rls_ simplify_Integral
   316 		   ],
   317 	 scr = EmptyScr};
   318 ruleset' := 
   319 overwritelthy thy (!ruleset', 
   320 	    [("integration_rules", prep_rls integration_rules),
   321 	     ("add_new_c", prep_rls add_new_c),
   322 	     ("simplify_Integral", prep_rls simplify_Integral),
   323 	     ("integration", prep_rls integration),
   324 	     ("separate_bdv2", separate_bdv2),
   325 	     ("norm_Rational_noadd_fractions", norm_Rational_noadd_fractions),
   326 	     ("norm_Rational_rls_noadd_fractions", 
   327 	      norm_Rational_rls_noadd_fractions)
   328 	     ]);
   329 
   330 (** problems **)
   331 
   332 store_pbt
   333  (prep_pbt (theory "Integrate") "pbl_fun_integ" [] e_pblID
   334  (["integrate","function"],
   335   [("#Given" ,["functionTerm f_", "integrateBy v_"]),
   336    ("#Find"  ,["antiDerivative F_"])
   337   ],
   338   append_rls "e_rls" e_rls [(*for preds in where_*)], 
   339   SOME "Integrate (f_, v_)", 
   340   [["diff","integration"]]));
   341  
   342 (*here "named" is used differently from Differentiation"*)
   343 store_pbt
   344  (prep_pbt (theory "Integrate") "pbl_fun_integ_nam" [] e_pblID
   345  (["named","integrate","function"],
   346   [("#Given" ,["functionTerm f_", "integrateBy v_"]),
   347    ("#Find"  ,["antiDerivativeName F_"])
   348   ],
   349   append_rls "e_rls" e_rls [(*for preds in where_*)], 
   350   SOME "Integrate (f_, v_)", 
   351   [["diff","integration","named"]]));
   352  
   353 (** methods **)
   354 
   355 store_met
   356     (prep_met (theory "Integrate") "met_diffint" [] e_metID
   357 	      (["diff","integration"],
   358 	       [("#Given" ,["functionTerm f_", "integrateBy v_"]),
   359 		("#Find"  ,["antiDerivative F_"])
   360 		],
   361 	       {rew_ord'="tless_true", rls'=Atools_erls, calc = [], 
   362 		srls = e_rls, 
   363 		prls=e_rls,
   364 	     crls = Atools_erls, nrls = e_rls},
   365 "Script IntegrationScript (f_::real) (v_::real) =                " ^
   366 "  (let t_ = Take (Integral f_ D v_)                             " ^
   367 "   in (Rewrite_Set_Inst [(bdv,v_)] integration False) (t_::real))"
   368 ));
   369     
   370 store_met
   371     (prep_met (theory "Integrate") "met_diffint_named" [] e_metID
   372 	      (["diff","integration","named"],
   373 	       [("#Given" ,["functionTerm f_", "integrateBy v_"]),
   374 		("#Find"  ,["antiDerivativeName F_"])
   375 		],
   376 	       {rew_ord'="tless_true", rls'=Atools_erls, calc = [], 
   377 		srls = e_rls, 
   378 		prls=e_rls,
   379 		crls = Atools_erls, nrls = e_rls},
   380 "Script NamedIntegrationScript (f_::real) (v_::real) (F_::real=>real) = " ^
   381 "  (let t_ = Take (F_ v_ = Integral f_ D v_)                            " ^
   382 "   in ((Try (Rewrite_Set_Inst [(bdv,v_)] simplify_Integral False)) @@  " ^
   383 "       (Rewrite_Set_Inst [(bdv,v_)] integration False)) t_)            "
   384     ));
   385 *}
   386 
   387 end