section ‹Auxiliary Lemmas›
theory ODE_Auxiliarities
imports
  "HOL-Analysis.Analysis"
  "HOL-Library.Float"
  "List-Index.List_Index"
  Affine_Arithmetic.Affine_Arithmetic_Auxiliarities
  Affine_Arithmetic.Executable_Euclidean_Space
begin

instantiation prod :: (zero_neq_one, zero_neq_one) zero_neq_one
begin

definition "1 = (1, 1)"

instance by standard (simp add: zero_prod_def one_prod_def)
end

subsection ‹there is no inner product for type @{typ "'a ⇒⇩L 'b"}›

lemma (in real_inner) parallelogram_law: "(norm (x + y))⇧2 + (norm (x - y))⇧2 = 2 * (norm x)⇧2 + 2 * (norm y)⇧2"
proof -
  have "(norm (x + y))⇧2 + (norm (x - y))⇧2 = inner (x + y) (x + y) + inner (x - y) (x - y)"
    by (simp add: norm_eq_sqrt_inner)
  also have "… = 2 * (norm x)⇧2 + 2 * (norm y)⇧2"
    by (simp add: algebra_simps norm_eq_sqrt_inner)
  finally show ?thesis .
qed

locale no_real_inner
begin

lift_definition fstzero::"(real*real) ⇒⇩L (real*real)" is "λ(x, y). (x, 0)"
  by (auto intro!: bounded_linearI')

lemma [simp]: "fstzero (a, b) = (a, 0)"
  by transfer simp

lift_definition zerosnd::"(real*real) ⇒⇩L (real*real)" is "λ(x, y). (0, y)"
  by (auto intro!: bounded_linearI')

lemma [simp]: "zerosnd (a, b) = (0, b)"
  by transfer simp

lemma fstzero_add_zerosnd: "fstzero + zerosnd = id_blinfun"
  by transfer auto

lemma norm_fstzero_zerosnd: "norm fstzero = 1" "norm zerosnd = 1" "norm (fstzero - zerosnd) = 1"
  by (rule norm_blinfun_eqI[where x="(1, 0)"]) (auto simp: norm_Pair blinfun.bilinear_simps
    intro: norm_blinfun_eqI[where x="(0, 1)"] norm_blinfun_eqI[where x="(1, 0)"])

text ‹compare with @{thm parallelogram_law}›

lemma "(norm (fstzero + zerosnd))⇧2 + (norm (fstzero - zerosnd))⇧2 ≠
    2 * (norm fstzero)⇧2 + 2 * (norm zerosnd)⇧2"
  by (simp add: fstzero_add_zerosnd norm_fstzero_zerosnd)

end

subsection ‹Topology›

subsection ‹Vector Spaces›

lemma ex_norm_eq_1: "∃x. norm (x::'a::{real_normed_vector, perfect_space}) = 1"
  by (metis vector_choose_size zero_le_one)

subsection ‹Reals›

subsection ‹Balls›

text ‹sometimes @{thm mem_ball} etc. are not good ‹[simp]› rules (although they are often useful):
  not sure that inequalities are ``simpler'' than set membership (distorts automatic reasoning
  when only sets are involved)›
lemmas [simp del] = mem_ball mem_cball mem_sphere mem_ball_0 mem_cball_0


subsection ‹Boundedness›

lemma bounded_subset_cboxE:
  assumes "⋀i. i ∈ Basis ⟹ bounded ((λx. x ∙ i) ` X)"
  obtains a b where "X ⊆ cbox a b"
proof -
  have "⋀i. i ∈ Basis ⟹ ∃a b. ((λx. x ∙ i) ` X) ⊆ {a..b}"
    by (metis box_real(2) box_subset_cbox subset_trans bounded_subset_box_symmetric[OF assms] )
  then obtain a b where bnds: "⋀i. i ∈ Basis ⟹ ((λx. x ∙ i) ` X) ⊆ {a i .. b i}" 
    by metis
  then have "X ⊆ {x. ∀i∈Basis. x ∙ i ∈ {a i .. b i}}"
    by force
  also have "… = cbox (∑i∈Basis. a i *⇩R i) (∑i∈Basis. b i *⇩R i)"
    by (auto simp: cbox_def)
  finally show ?thesis ..
qed

lemma
  bounded_euclideanI:
  assumes "⋀i. i ∈ Basis ⟹ bounded ((λx. x ∙ i) ` X)"
  shows "bounded X"
proof -
  from bounded_subset_cboxE[OF assms] obtain a b where "X ⊆ cbox a b" .
  with bounded_cbox show ?thesis by (rule bounded_subset)
qed

subsection ‹Intervals›

notation closed_segment ("(1{_--_})")
notation open_segment ("(1{_<--<_})")

lemma min_zero_mult_nonneg_le: "0 ≤ h' ⟹ h' ≤ h ⟹ min 0 (h * k::real) ≤ h' * k"
  by (metis dual_order.antisym le_cases min_le_iff_disj mult_eq_0_iff mult_le_0_iff
      mult_right_mono_neg)

lemma max_zero_mult_nonneg_le: "0 ≤ h' ⟹ h' ≤ h ⟹ h' * k ≤ max 0 (h * k::real)"
  by (metis dual_order.antisym le_cases le_max_iff_disj mult_eq_0_iff mult_right_mono
      zero_le_mult_iff)

lemmas closed_segment_eq_real_ivl = closed_segment_eq_real_ivl

lemma bdd_above_is_intervalI: "bdd_above I" if "is_interval I" "a ≤ b" "a ∈ I" "b ∉ I" for I::"real set"
  by (meson bdd_above_def is_interval_1 le_cases that)

lemma bdd_below_is_intervalI: "bdd_below I" if "is_interval I" "a ≤ b" "a ∉ I" "b ∈ I" for I::"real set"
  by (meson bdd_below_def is_interval_1 le_cases that)


subsection ‹Extended Real Intervals›

subsection ‹Euclidean Components›

subsection ‹Operator Norm›

subsection ‹Limits›

lemma eventually_open_cball:
  assumes "open X"
  assumes "x ∈ X"
  shows "eventually (λe. cball x e ⊆ X) (at_right 0)"
proof -
  from open_contains_cball_eq[OF assms(1)] assms(2)
  obtain e where "e > 0" "cball x e ⊆ X" by auto
  thus ?thesis
    by (auto simp: eventually_at dist_real_def mem_cball intro!: exI[where x=e])
qed

subsection ‹Continuity›

subsection ‹Derivatives›

lemma
  if_eventually_has_derivative:
  assumes "(f has_derivative F') (at x within S)"
  assumes "∀⇩F x in at x within S. P x" "P x" "x ∈ S"
  shows "((λx. if P x then f x else g x) has_derivative F') (at x within S)"
  using assms(1)
  apply (rule has_derivative_transform_eventually)
  subgoal using assms(2) by eventually_elim auto
  by (auto simp: assms)


lemma norm_le_in_cubeI: "norm x ≤ norm y"
  if "⋀i. i ∈ Basis ⟹ abs (x ∙ i) ≤ abs (y ∙ i)" for x y
  unfolding norm_eq_sqrt_inner
  apply (subst euclidean_inner)
  apply (subst (3) euclidean_inner)
  using that
  by (auto intro!: sum_mono simp: abs_le_square_iff power2_eq_square[symmetric])

lemma has_derivative_partials_euclidean_convexI:
  fixes f::"'a::euclidean_space ⇒ 'b::real_normed_vector"
  assumes f': "⋀i x xi. i ∈ Basis ⟹ (∀j∈Basis. x ∙ j ∈ X j) ⟹ xi = x ∙ i ⟹
    ((λp. f (x + (p - x ∙ i) *⇩R i)) has_vector_derivative f' i x) (at xi within X i)"
  assumes df_cont: "⋀i. i ∈ Basis ⟹ (f' i ⤏ (f' i x)) (at x within {x. ∀j∈Basis. x ∙ j ∈ X j})"
  assumes "⋀i. i ∈ Basis ⟹ x ∙ i ∈ X i"
  assumes "⋀i. i ∈ Basis ⟹ convex (X i)"
  shows "(f has_derivative (λh. ∑j∈Basis. (h ∙ j) *⇩R f' j x)) (at x within {x. ∀j∈Basis. x ∙ j ∈ X j})"
    (is "_ (at x within ?S)")
proof (rule has_derivativeI)
  show "bounded_linear (λh. ∑j∈Basis. (h ∙ j) *⇩R f' j x)"
    by (auto intro!: bounded_linear_intros)

  obtain E where [simp]: "set E = (Basis::'a set)" "distinct E"
    using finite_distinct_list[OF finite_Basis] by blast

  have [simp]: "length E = DIM('a)"
    using ‹distinct E› distinct_card by fastforce
  have [simp]: "E ! j ∈ Basis" if "j < DIM('a)" for j
    by (metis ‹length E = DIM('a)› ‹set E = Basis› nth_mem that)
  have "convex ?S"
    by (rule convex_prod) (use assms in auto)

  have sum_Basis_E: "sum g Basis = (∑j<DIM('a). g (E ! j))" for g
    apply (rule sum.reindex_cong[OF _ _ refl])
    apply (auto simp: inj_on_nth)
    by (metis ‹distinct E› ‹length E = DIM('a)› ‹set E = Basis› bij_betw_def bij_betw_nth)

  have segment: "∀⇩F x' in at x within ?S. x' ∈ ?S" "∀⇩F x' in at x within ?S. x' ≠ x"
    unfolding eventually_at_filter by auto


  show "((λy. (f y - f x - (∑j∈Basis. ((y - x) ∙ j) *⇩R f' j x)) /⇩R norm (y - x)) ⤏ 0) (at x within {x. ∀j∈Basis. x ∙ j ∈ X j})"
    apply (rule tendstoI)
    unfolding norm_conv_dist[symmetric]
  proof -
    fix e::real
    assume "e > 0"
    define B where "B = e / norm (2*DIM('a) + 1)"
    with ‹e > 0› have B_thms: "B > 0" "2 * DIM('a) * B < e" "B ≥ 0"
      by (auto simp: divide_simps algebra_simps B_def)
    define B' where "B' = B / 2"
    have "B' > 0" by (simp add: B'_def ‹0 < B›)
    have "∀i ∈ Basis. ∀⇩F xa in at x within {x. ∀j∈Basis. x ∙ j ∈ X j}. dist (f' i xa) (f' i x) < B'"
      apply (rule ballI)
      subgoal premises prems using df_cont[OF prems, THEN tendstoD, OF ‹0 < B'›] .
      done
    from eventually_ball_finite[OF finite_Basis this]
    have "∀⇩F x' in at x within {x. ∀j∈Basis. x ∙ j ∈ X j}. ∀j∈Basis. dist (f' j x') (f' j x) < B'" .
    then obtain d where "d > 0"
      and "⋀x' j. x' ∈ {x. ∀j∈Basis. x ∙ j ∈ X j} ⟹ x' ≠ x ⟹ dist x' x < d ⟹ j ∈ Basis ⟹ dist (f' j x') (f' j x) < B'"
      using ‹0 < B'›
      by (auto simp: eventually_at)
    then have B': "x' ∈ {x. ∀j∈Basis. x ∙ j ∈ X j} ⟹ dist x' x < d ⟹ j ∈ Basis ⟹ dist (f' j x') (f' j x) < B'" for x' j
      by (cases "x' = x", auto simp add: ‹0 < B'›)
    then have B: "norm (f' j x' - f' j y) < B" if
      "(⋀j. j ∈ Basis ⟹ x' ∙ j ∈ X j)"
      "(⋀j. j ∈ Basis ⟹ y ∙ j ∈ X j)"
      "dist x' x < d"
      "dist y x < d"
      "j ∈ Basis"
      for x' y j
    proof -
      have "dist (f' j x') (f' j x) < B'" "dist (f' j y) (f' j x) < B'"
        using that
        by (auto intro!: B')
      then have "dist (f' j x') (f' j y) < B' + B'" by norm
      also have "… = B" by (simp add: B'_def)
      finally show ?thesis by (simp add: dist_norm)
    qed
    have "∀⇩F x' in at x within {x. ∀j∈Basis. x ∙ j ∈ X j}. dist x' x < d"
      by (rule tendstoD[OF tendsto_ident_at ‹d > 0›])
    with segment
    show "∀⇩F x' in at x within {x. ∀j∈Basis. x ∙ j ∈ X j}.
      norm ((f x' - f x - (∑j∈Basis. ((x' - x) ∙ j) *⇩R f' j x)) /⇩R norm (x' - x)) < e"
    proof eventually_elim
      case (elim x')
      then have os_subset: "open_segment x x' ⊆ ?S"
        using ‹convex ?S› assms(3)
        unfolding convex_contains_open_segment
        by auto
      then have cs_subset: "closed_segment x x' ⊆ ?S"
        using elim assms(3) by (auto simp add: open_segment_def)
      have csc_subset: "closed_segment (x' ∙ i) (x ∙ i) ⊆ X i" if i: "i ∈ Basis" for i
        apply (rule closed_segment_subset)
        using cs_subset elim assms(3,4) that
        by (auto )
      have osc_subset: "open_segment (x' ∙ i) (x ∙ i) ⊆ X i" if i: "i ∈ Basis" for i
        using segment_open_subset_closed csc_subset[OF i]
        by (rule order_trans)

      define h where "h = x' - x"
      define z where "z j = (∑k<j. (h ∙ E ! k) *⇩R (E ! k))" for j
      define g where "g j t = (f (x + z j + (t - x ∙ E ! j) *⇩R E ! j))" for j t
      have z: "z j ∙ E ! j = 0" if "j < DIM('a)" for j
        using that
        by (auto simp: z_def h_def algebra_simps inner_sum_left inner_Basis if_distrib
            nth_eq_iff_index_eq
            sum.delta
            intro!: euclidean_eqI[where 'a='a]
            cong: if_cong)
      from distinct_Ex1[OF ‹distinct E›, unfolded ‹set E = Basis› Ex1_def ‹length E = _›]
      obtain index where
        index: "⋀i. i ∈ Basis ⟹ i = E ! index i" "⋀i. i ∈ Basis ⟹ index i < DIM('a)"
        and unique: "⋀i j. i ∈ Basis ⟹ j < DIM('a) ⟹ E ! j = i ⟹ j = index i"
        by metis
      have nth_eq_iff_index: "E ! k = i ⟷ index i = k" if "i ∈ Basis" "k < DIM('a)" for k i
        using unique[OF that] index[OF ‹i ∈ Basis›]
        by auto
      have z_inner: "z j ∙ i = (if j ≤ index i then 0 else h ∙ i)" if "j < DIM('a)" "i ∈ Basis" for j i
        using that index[of i]
        by (auto simp: z_def h_def algebra_simps inner_sum_left inner_Basis if_distrib
            sum.delta nth_eq_iff_index
            intro!: euclidean_eqI[where 'a='a]
            cong: if_cong)
      have z_mem: "j < DIM('a) ⟹ ja ∈ Basis ⟹ x ∙ ja + z j ∙ ja ∈ X ja" for j ja
        using csc_subset
        by (auto simp: z_inner h_def algebra_simps)
      have "norm (x' - x) < d"
        using elim by (simp add: dist_norm)
      have norm_z': "y ∈ closed_segment (x ∙ E ! j) (x' ∙ E ! j) ⟹ norm (z j + y *⇩R E ! j - (x ∙ E ! j) *⇩R E ! j) < d"
        if "j < DIM('a)"
        for j y
        apply (rule le_less_trans[OF _ ‹norm (x' - x) < d›])
        apply (rule norm_le_in_cubeI)
        apply (auto simp: inner_diff_left inner_add_left inner_Basis that z)
        subgoal by (auto simp: closed_segment_eq_real_ivl split: if_splits)
        subgoal for i
          using that
          by (auto simp: z_inner h_def algebra_simps)
        done
      have norm_z: "norm (z j) < d" if "j < DIM('a)" for j
        using norm_z'[OF that ends_in_segment(1)]
        by (auto simp: z_def)
      {
        fix j
        assume j: "j < DIM('a)"
        have eq: "(x + z j + ((y - (x + z j)) ∙ E ! j) *⇩R E ! j +
          (p - (x + z j + ((y - (x + z j)) ∙ E ! j) *⇩R E ! j) ∙ E ! j) *⇩R
          E ! j) = (x + z j + (p - (x ∙ E ! j)) *⇩R E ! j)" for y p
          by (auto simp: algebra_simps j z)
        have f_has_derivative: "((λp. f (x + z j + (p - x ∙ E ! j) *⇩R E ! j)) has_derivative (λxa. xa *⇩R f' (E ! j) (x + z j + ((y *⇩R E ! j - (x + z j)) ∙ E ! j) *⇩R E ! j)))
            (at y within closed_segment (x ∙ E ! j) (x' ∙ E ! j))"
          if "y ∈ closed_segment (x ∙ E ! j) (x' ∙ E ! j)"
          for y
          apply (rule has_derivative_subset)
           apply (rule f'[unfolded has_vector_derivative_def,
                where x= "x + z j + ((y *⇩R E!j - (x + z j))∙ E!j) *⇩R E ! j" and i="E ! j", unfolded eq])
          subgoal by (simp add: j)
          subgoal
            using that
            apply (auto simp: algebra_simps z j inner_Basis)
            using closed_segment_commute ‹E ! j ∈ Basis› csc_subset apply blast
            by (simp add: z_mem j)
          subgoal by (auto simp: algebra_simps z j inner_Basis)
          subgoal
            apply (auto simp: algebra_simps z j inner_Basis)
            using closed_segment_commute ‹⋀j. j < DIM('a) ⟹ E ! j ∈ Basis› csc_subset j apply blast
            done
          done
        have *: "((xa *⇩R E ! j - (x + z j)) ∙ E ! j) = xa - x ∙ E ! j" for xa
          by (auto simp: algebra_simps z j)
        have g': "(g j has_vector_derivative (f' (E ! j) (x + z j + (xa - x ∙ E ! j) *⇩R E ! j)))
          (at xa within (closed_segment (x∙E!j) (x'∙E!j)))"
          (is "(_ has_vector_derivative ?g' j xa) _")
          if "xa ∈ closed_segment (x∙E!j) (x'∙E!j)" for xa
          using that
          by (auto simp: has_vector_derivative_def g_def[abs_def] *
              intro!: derivative_eq_intros f_has_derivative[THEN has_derivative_eq_rhs])
        define g' where "g' j = ?g' j" for j
        with g' have g': "(g j has_vector_derivative g' j t) (at t within closed_segment (x∙E!j) (x'∙E!j))"
          if "t ∈ closed_segment (x∙E!j) (x'∙E!j)"
          for t
          by (simp add: that)
        have cont_bound: "⋀y. y∈closed_segment (x ∙ E ! j) (x' ∙ E ! j) ⟹ norm (g' j y - g' j (x ∙ E ! j)) ≤ B"
          apply (auto simp add: g'_def j algebra_simps inner_Basis z dist_norm
              intro!: less_imp_le B z_mem norm_z norm_z')
          using closed_segment_commute ‹⋀j. j < DIM('a) ⟹ E ! j ∈ Basis› csc_subset j apply blast
          done
        from vector_differentiable_bound_linearization[OF g' order_refl cont_bound ends_in_segment(1)]
        have n: "norm (g j (x' ∙ E ! j) - g j (x ∙ E ! j) - (x' ∙ E ! j - x ∙ E ! j) *⇩R g' j (x ∙ E ! j)) ≤ norm (x' ∙ E ! j - x ∙ E ! j) * B"
          .
        have "z (Suc j) = z j + (x' ∙ E ! j - x ∙ E ! j) *⇩R E ! j"
          by (auto simp: z_def h_def algebra_simps)
        then have "f (x + z (Suc j)) = f (x + z j + (x' ∙ E ! j - x ∙ E ! j) *⇩R E ! j) "
          by (simp add: ac_simps)
        with n have "norm (f (x + z (Suc j)) - f (x + z j) - (x' ∙ E ! j - x ∙ E ! j) *⇩R f' (E ! j) (x + z j)) ≤ ¦x' ∙ E ! j - x ∙ E ! j¦ * B"
          by (simp add: g_def g'_def)
      } note B_le = this
      have B': "norm (f' (E ! j) (x + z j) - f' (E ! j) x) ≤ B" if "j < DIM('a)" for j
        using that assms(3)
        by (auto simp add: algebra_simps inner_Basis z dist_norm ‹0 < d›
            intro!: less_imp_le B z_mem norm_z)
      have "(∑j≤DIM('a) - 1. f (x + z j) - f (x + z (Suc j))) = f (x + z 0) - f (x + z (Suc (DIM('a) - 1)))"
        by (rule sum_telescope)
      moreover have "z DIM('a) = h"
        using index
        by (auto simp: z_def h_def algebra_simps inner_sum_left inner_Basis if_distrib
            nth_eq_iff_index 
            sum.delta 
            intro!: euclidean_eqI[where 'a='a]
            cong: if_cong)
      moreover have "z 0 = 0"
        by (auto simp: z_def)
      moreover have "{..DIM('a) - 1} = {..<DIM('a)}"
        using le_imp_less_Suc by fastforce
      ultimately have "f x - f (x + h) = (∑j<DIM('a). f (x + z j) - f (x + z (Suc j)))"
        by (auto simp: )
      then have "norm (f (x + h) - f x - (∑j∈Basis. ((x' - x) ∙ j) *⇩R f' j x)) =
        norm(
          (∑j<DIM('a). f (x + z (Suc j)) - f (x + z j) - (x' ∙ E ! j - x ∙ E ! j) *⇩R f' (E ! j) (x + z j)) +
          (∑j<DIM('a). (x' ∙ E ! j - x ∙ E ! j) *⇩R (f' (E ! j) (x + z j) - f' (E ! j) x)))"
        (is "_ = norm (sum ?a ?E + sum ?b ?E)")
        by (intro arg_cong[where f=norm]) (simp add: sum_negf sum_subtractf sum.distrib algebra_simps sum_Basis_E)
      also have "… ≤ norm (sum ?a ?E) + norm (sum ?b ?E)" by (norm)
      also have "norm (sum ?a ?E) ≤ sum (λx. norm (?a x)) ?E"
        by (rule norm_sum)
      also have "… ≤ sum (λj. norm ¦x' ∙ E ! j - x ∙ E ! j¦ * B) ?E"
        by (auto intro!: sum_mono B_le)
      also have "… ≤ sum (λj. norm (x' - x) * B) ?E"
        apply (rule sum_mono)
        apply (auto intro!: mult_right_mono ‹0 ≤ B›)
        by (metis (full_types) ‹⋀j. j < DIM('a) ⟹ E ! j ∈ Basis› inner_diff_left norm_bound_Basis_le order_refl)
      also have "… = norm (x' - x) * DIM('a) * B"
        by simp
      also have "norm (sum ?b ?E) ≤ sum (λx. norm (?b x)) ?E"
        by (rule norm_sum)
      also have "… ≤ sum (λj. norm (x' - x) * B) ?E"
        apply (intro sum_mono)
        apply (auto intro!: mult_mono B')
         apply (metis (full_types) ‹⋀j. j < DIM('a) ⟹ E ! j ∈ Basis› inner_diff_left norm_bound_Basis_le order_refl)
        done
      also have "… = norm (x' - x) * DIM('a) * B"
        by simp
      finally have "norm (f (x + h) - f x - (∑j∈Basis. ((x' - x) ∙ j) *⇩R f' j x)) ≤
          norm (x' - x) * real DIM('a) * B + norm (x' - x) * real DIM('a) * B"
        by arith
      also have "… /⇩R norm (x' - x) ≤ norm (2 * DIM('a) * B)"
        using ‹B ≥ 0›
        by (simp add: divide_simps abs_mult)
      also have "… < e" using B_thms by simp
      finally show ?case by (auto simp: divide_simps abs_mult h_def)
    qed
  qed
qed

lemma
  frechet_derivative_equals_partial_derivative:
  fixes f::"'a::euclidean_space ⇒ 'a"
  assumes Df: "⋀x. (f has_derivative Df x) (at x)"
  assumes f': "((λp. f (x + (p - x ∙ i) *⇩R i) ∙ b) has_real_derivative f' x i b) (at (x ∙ i))"
  shows "Df x i ∙ b = f' x i b"
proof -
  define Dfb where "Dfb x = Blinfun (Df x)" for x
  have Dfb_apply: "blinfun_apply (Dfb x) = Df x" for x
    unfolding Dfb_def
    apply (rule bounded_linear_Blinfun_apply)
    apply (rule has_derivative_bounded_linear)
    apply (rule assms)
    done
  have "Dfb x = blinfun_of_matrix (λi b. Dfb x b ∙ i)" for x
    using blinfun_of_matrix_works[of "Dfb x"] by auto
  have Dfb: "⋀x. (f has_derivative Dfb x) (at x)"
    by (auto simp: Dfb_apply Df)
  note [derivative_intros] = diff_chain_at[OF _ Dfb, unfolded o_def]
  have "((λp. f (x + (p - x ∙ i) *⇩R i) ∙ b) has_real_derivative Dfb x i ∙ b) (at (x ∙ i))"
    by (auto intro!: derivative_eq_intros ext simp: has_field_derivative_def blinfun.bilinear_simps)
  from DERIV_unique[OF f' this]
  show ?thesis by (simp add: Dfb_apply)
qed


subsection ‹Integration›

lemmas content_real[simp]
lemmas integrable_continuous[intro, simp]
  and integrable_continuous_real[intro, simp]


lemma integral_eucl_le:
  fixes f g::"'a::euclidean_space ⇒ 'b::ordered_euclidean_space"
  assumes "f integrable_on s"
    and "g integrable_on s"
    and "⋀x. x ∈ s ⟹ f x ≤ g x"
  shows "integral s f ≤ integral s g"
  using assms
  by (auto intro!: integral_component_le simp: eucl_le[where 'a='b])

lemma
  integral_ivl_bound:
  fixes l u::"'a::ordered_euclidean_space"
  assumes "⋀x h'. h' ∈ {t0 .. h} ⟹ x ∈ {t0 .. h} ⟹ (h' - t0) *⇩R f x ∈ {l .. u}"
  assumes "t0 ≤ h"
  assumes f_int: "f integrable_on {t0 .. h}"
  shows "integral {t0 .. h} f ∈ {l .. u}"
proof -
  from assms(1)[of t0 t0] assms(2) have "0 ∈ {l .. u}" by auto
  have "integral {t0 .. h} f = integral {t0 .. h} (λt. if t ∈ {t0, h} then 0 else f t)"
    by (rule integral_spike[where S="{t0, h}"]) auto
  also
  {
    fix x
    assume 1: "x ∈ {t0 <..< h}"
    with assms have "(h - t0) *⇩R f x ∈ {l .. u}" by auto
    then have "(if x ∈ {t0, h} then 0 else f x) ∈ {l /⇩R (h - t0) .. u /⇩R (h - t0)}"
      using ‹x ∈ _›
      by (auto simp: inverse_eq_divide
        simp: eucl_le[where 'a='a] field_simps algebra_simps)
  }
  then have "… ∈ {integral {t0..h} (λ_. l /⇩R (h - t0)) .. integral {t0..h} (λ_. u /⇩R (h - t0))}"
    unfolding atLeastAtMost_iff
    apply (safe intro!: integral_eucl_le)
    using ‹0 ∈ {l .. u}›
    apply (auto intro!: assms
      intro: integrable_continuous_real  integrable_spike[where S="{t0, h}", OF f_int]
      simp: eucl_le[where 'a='a] divide_simps
      split: if_split_asm)
    done
  also have "… ⊆ {l .. u}"
    using assms ‹0 ∈ {l .. u}›
    by (auto simp: inverse_eq_divide)
  finally show ?thesis .
qed

lemma
  add_integral_ivl_bound:
  fixes l u::"'a::ordered_euclidean_space"
  assumes "⋀x h'. h' ∈ {t0 .. h} ⟹ x ∈ {t0 .. h} ⟹ (h' - t0) *⇩R f x ∈ {l - x0 .. u - x0}"
  assumes "t0 ≤ h"
  assumes f_int: "f integrable_on {t0 .. h}"
  shows "x0 + integral {t0 .. h} f ∈ {l .. u}"
  using integral_ivl_bound[OF assms]
  by (auto simp: algebra_simps)

subsection ‹conditionally complete lattice›


subsection ‹Lists›

lemma
  Ball_set_Cons[simp]: "(∀a∈set_Cons x y. P a) ⟷ (∀a∈x. ∀b∈y. P (a#b))"
  by (auto simp: set_Cons_def)

lemma set_cons_eq_empty[iff]: "set_Cons a b = {} ⟷ a = {} ∨ b = {}"
  by (auto simp: set_Cons_def)

lemma listset_eq_empty_iff[iff]: "listset XS = {} ⟷ {} ∈ set XS"
  by (induction XS) auto

lemma sing_in_sings[simp]: "[x] ∈ (λx. [x]) ` xd ⟷ x ∈ xd"
  by auto

lemma those_eq_None_set_iff: "those xs = None ⟷ None ∈ set xs"
  by (induction xs) (auto split: option.split)

lemma those_eq_Some_lengthD: "those xs = Some ys ⟹ length xs = length ys"
  by (induction xs arbitrary: ys) (auto split: option.splits)

lemma those_eq_Some_map_Some_iff: "those xs = Some ys ⟷ (xs = map Some ys)" (is "?l ⟷ ?r")
proof safe
  assume ?l
  then have "length xs = length ys"
    by (rule those_eq_Some_lengthD)
  then show ?r using ‹?l›
    by (induction xs ys rule: list_induct2) (auto split: option.splits)
next
  assume ?r
  then have "length xs = length ys"
    by simp
  then show "those (map Some ys) = Some ys" using ‹?r›
    by (induction xs ys rule: list_induct2) (auto split: option.splits)
qed


subsection ‹Set(sum)›

subsection ‹Max›

subsection ‹Uniform Limit›

subsection ‹Bounded Linear Functions›

lift_definition comp3::― ‹TODO: name?›
  "('c::real_normed_vector ⇒⇩L 'd::real_normed_vector) ⇒ ('b::real_normed_vector ⇒⇩L 'c) ⇒⇩L 'b ⇒⇩L 'd" is
  "λ(cd::('c ⇒⇩L 'd)) (bc::'b ⇒⇩L 'c). (cd o⇩L bc)"
  by (rule bounded_bilinear.bounded_linear_right[OF bounded_bilinear_blinfun_compose])

lemma blinfun_apply_comp3[simp]: "blinfun_apply (comp3 a) b = (a o⇩L b)"
  by (simp add: comp3.rep_eq)

lemma bounded_linear_comp3[bounded_linear]: "bounded_linear comp3"
  by transfer (rule bounded_bilinear_blinfun_compose)

lift_definition comp12::― ‹TODO: name?›
  "('a::real_normed_vector ⇒⇩L 'c::real_normed_vector) ⇒ ('b::real_normed_vector ⇒⇩L 'c) ⇒ ('a × 'b) ⇒⇩L 'c"
  is "λf g (a, b). f a + g b"
  by (auto intro!: bounded_linear_intros
    intro: bounded_linear_compose
    simp: split_beta')

lemma blinfun_apply_comp12[simp]: "blinfun_apply (comp12 f g) b = f (fst b) + g (snd b)"
  by (simp add: comp12.rep_eq split_beta)


subsection ‹Order Transitivity Attributes›

attribute_setup le = ‹Scan.succeed (Thm.rule_attribute [] (fn context => fn thm => thm RS @{thm order_trans}))›
  "transitive version of inequality (useful for intro)"
attribute_setup ge = ‹Scan.succeed (Thm.rule_attribute [] (fn context => fn thm => thm RS @{thm order_trans[rotated]}))›
  "transitive version of inequality (useful for intro)"


subsection ‹point reflection›

definition preflect::"'a::real_vector ⇒ 'a ⇒ 'a" where "preflect ≡ λt0 t. 2 *⇩R t0 - t"

lemma preflect_preflect[simp]: "preflect t0 (preflect t0 t) = t"
  by (simp add: preflect_def algebra_simps)

lemma preflect_preflect_image[simp]: "preflect t0 ` preflect t0 ` S = S"
  by (simp add: image_image)

lemma is_interval_preflect[simp]: "is_interval (preflect t0 ` S) ⟷ is_interval S"
  by (auto simp: preflect_def[abs_def])

lemma iv_in_preflect_image[intro, simp]: "t0 ∈ T ⟹ t0 ∈ preflect t0 ` T"
  by (auto intro!: image_eqI simp: preflect_def algebra_simps scaleR_2)

lemma preflect_tendsto[tendsto_intros]:
  fixes l::"'a::real_normed_vector"
  shows "(g ⤏ l) F ⟹ (h ⤏ m) F ⟹ ((λx. preflect (g x) (h x)) ⤏ preflect l m) F"
  by (auto intro!: tendsto_eq_intros simp: preflect_def)

lemma continuous_preflect[continuous_intros]:
  fixes a::"'a::real_normed_vector"
  shows "continuous (at a within A) (preflect t0)"
  by (auto simp: continuous_within intro!: tendsto_intros)

lemma
  fixes t0::"'a::ordered_real_vector"
  shows preflect_le[simp]: "t0 ≤ preflect t0 b ⟷ b ≤ t0"
    and le_preflect[simp]: "preflect t0 b ≤ t0 ⟷ t0 ≤ b"
    and antimono_preflect: "antimono (preflect t0)"
    and preflect_le_preflect[simp]: "preflect t0 a ≤ preflect t0 b ⟷ b ≤ a"
    and preflect_eq_cancel[simp]: "preflect t0 a = preflect t0 b ⟷ a = b"
  by (auto intro!: antimonoI simp: preflect_def scaleR_2)

lemma preflect_eq_point_iff[simp]: "t0 = preflect t0 s ⟷ t0 = s" "preflect t0 s = t0 ⟷ t0 = s"
  by (auto simp: preflect_def algebra_simps scaleR_2)

lemma preflect_minus_self[simp]: "preflect t0 s - t0 = t0 - s"
  by (simp add: preflect_def scaleR_2)

end

theory MVT_Ex
imports
  "HOL-Analysis.Analysis"
  "HOL-Decision_Procs.Approximation"
  "../ODE_Auxiliarities"
begin

subsection ‹(Counter)Example of Mean Value Theorem in Euclidean Space \label{sec:countermvt}›

text ‹There is no exact analogon of the mean value theorem in the multivariate case!›

lemma MVT_wrong: assumes
  "⋀J a u (f::real*real⇒real*real).
      (⋀x. FDERIV f x :> J x) ⟹
      (∃t∈{0<..<1}. f (a + u) - f a = J (a + t *⇩R u) u)"
  shows "False"
proof -
  have "⋀t::real*real. FDERIV (λt. (cos (fst t), sin (fst t))) t :> (λh. (- ((fst h) * sin (fst t)), (fst h) * cos (fst t)))"
    by (auto intro!: derivative_eq_intros)
  from assms[OF this, of "(pi, pi)" "(pi, pi)"] obtain t::real where t: "0 < t" "t < 1" and
    "pi * sin (t * pi) = 2" "cos (t * pi) = 0"
    by auto
  then obtain n where tpi: "t * pi = real_of_int n * (pi / 2)" and "odd n"
    by (auto simp: cos_zero_iff_int)
  then have teq: "t = real_of_int n / 2" by auto
  then have "n = 1" using t ‹odd n› by arith
  then have "t = 1/2" using teq by simp
  have "sin (t * pi) = 1"
    by (simp add: ‹t = 1/2› sin_eq_1)
  with ‹pi * sin (t * pi) = 2›
  have "pi = 2" by simp
  moreover have "pi > 2" using pi_approx by simp
  ultimately show False by simp
qed

lemma MVT_corrected:
  fixes f::"'a::ordered_euclidean_space⇒'b::euclidean_space"
  assumes fderiv: "⋀x. x ∈ D ⟹ (f has_derivative J x) (at x within D)"
  assumes line_in: "⋀x. ⟦0 ≤ x; x ≤ 1⟧ ⟹ a + x *⇩R u ∈ D"
  shows "(∃t∈Basis→{0<..<1}. (f (a + u) - f a) = (∑i∈Basis. (J (a + t i *⇩R u) u ∙ i) *⇩R i))"
proof -
  {
    fix i::'b
    assume "i ∈ Basis"
    have subset: "((λx. a + x *⇩R u) ` {0..1}) ⊆ D"
      using line_in by force
    have "⋀x. ⟦0 ≤ x; x ≤ 1⟧ ⟹ ((λb. f (a + b *⇩R u) ∙ i) has_derivative (λb. b *⇩R J (a + x *⇩R u) u ∙ i)) (at x within {0..1})"
      using line_in
      by (auto intro!: derivative_eq_intros
        has_derivative_subset[OF _ subset]
        has_derivative_in_compose[where f="λx. a + x *⇩R u"]
        fderiv line_in
        simp add: linear.scaleR[OF has_derivative_linear[OF fderiv]])
    with zero_less_one
    have "∃x∈{0<..<1}. f (a + 1 *⇩R u) ∙ i - f (a + 0 *⇩R u) ∙ i = (1 - 0) *⇩R J (a + x *⇩R u) u ∙ i"
      by (rule mvt_simple)
  }
  then obtain t where "∀i∈Basis. t i ∈ {0<..<1} ∧ f (a + u) ∙ i - f a ∙ i = J (a + t i *⇩R u) u ∙ i"
    by atomize_elim (force intro!: bchoice)
  hence "t ∈ Basis → {0<..<1}" "⋀i. i ∈ Basis ⟹ (f (a + u) - f a) ∙ i = J (a + t i *⇩R u) u ∙ i"
    by (auto simp: inner_diff_left)
  moreover hence "(f (a + u) - f a) = (∑i∈Basis. (J (a + t i *⇩R u) u ∙ i) *⇩R i)"
    by (intro euclidean_eqI[where 'a='b]) simp
  ultimately show ?thesis by blast
qed

lemma MVT_ivl:
  fixes f::"'a::ordered_euclidean_space⇒'b::ordered_euclidean_space"
  assumes fderiv: "⋀x. x ∈ D ⟹ (f has_derivative J x) (at x within D)"
  assumes J_ivl: "⋀x. x ∈ D ⟹ J x u ∈ {J0 .. J1}"
  assumes line_in: "⋀x. x ∈ {0..1} ⟹ a + x *⇩R u ∈ D"
  shows "f (a + u) - f a ∈ {J0..J1}"
proof -
  from MVT_corrected[OF fderiv line_in] obtain t where
    t: "t∈Basis → {0<..<1}" and
    mvt: "f (a + u) - f a = (∑i∈Basis. (J (a + t i *⇩R u) u ∙ i) *⇩R i)"
    by auto
  note mvt
  also have "… ∈ {J0 .. J1}"
  proof -
    have J: "⋀i. i ∈ Basis ⟹ J0 ≤ J (a + t i *⇩R u) u"
            "⋀i. i ∈ Basis ⟹ J (a + t i *⇩R u) u ≤ J1"
      using J_ivl t line_in by (auto simp: Pi_iff)
    show ?thesis
      using J
      unfolding atLeastAtMost_iff eucl_le[where 'a='b]
      by auto
  qed
  finally show ?thesis .
qed

lemma MVT:
  shows
  "⋀J J0 J1 a u (f::real*real⇒real*real).
    (⋀x. FDERIV f x :> J x) ⟹
    (⋀x. J x u ∈ {J0 .. J1}) ⟹
    f (a + u) - f a ∈ {J0 .. J1}"
  by (rule_tac J = J in MVT_ivl[where D=UNIV]) auto

lemma MVT_ivl':
  fixes f::"'a::ordered_euclidean_space⇒'b::ordered_euclidean_space"
  assumes fderiv: "(⋀x. x ∈ D ⟹ (f has_derivative J x) (at x within D))"
  assumes J_ivl: "⋀x. x ∈ D ⟹ J x (a - b) ∈ {J0..J1}"
  assumes line_in: "⋀x. x ∈ {0..1} ⟹ b + x *⇩R (a - b) ∈ D"
  shows "f a ∈ {f b + J0..f b + J1}"
proof -
  have "f (b + (a - b)) - f b ∈ {J0 .. J1}"
    using J_ivl MVT_ivl fderiv line_in by blast
  thus ?thesis
    by (auto simp: diff_le_eq le_diff_eq ac_simps)
qed

end