Theory Interval_Utilities

(***********************************************************************************
 * Copyright (c) University of Exeter, UK
 *
 * All rights reserved.
 *
 * Redistribution and use in source and binary forms, with or without
 * modification, are permitted provided that the following conditions are met:
 *
 * * Redistributions of source code must retain the above copyright notice, this
 *
 * * Redistributions in binary form must reproduce the above copyright notice,
 *   this list of conditions and the following disclaimer in the documentation
 *   and/or other materials provided with the distribution.
 *
 * THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS "AS IS"
 * AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
 * IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE ARE
 * DISCLAIMED. IN NO EVENT SHALL THE COPYRIGHT HOLDER OR CONTRIBUTORS BE LIABLE
 * FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL
 * DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR
 * SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER
 * CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY,
 * OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE
 * OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
 *
 * SPDX-License-Identifier: BSD-3-Clause
 ***********************************************************************************)

chapter‹Interval Utilities›
theory
  Interval_Utilities
imports
  "HOL-Library.Interval"
  "HOL-Analysis.Analysis"
  "HOL-Library.Interval_Float"
begin

section‹Preliminaries›

lemma compact_set_of:
  fixes X::'a::{preorder,topological_space, ordered_euclidean_space} interval
  shows compact (set_of X)
  by(simp add:set_of_eq compact_interval[of "lower X" "upper X"])

lemma bounded_set_of:
  fixes X::'a::{preorder,topological_space, ordered_euclidean_space} interval
  shows bounded (set_of X)
  by(simp add:set_of_eq compact_interval[of "lower X" "upper X"])

lemma compact_img_set_of:
  fixes X :: real interval and f::real  real
  assumes continuous_on (set_of X) f
  shows compact (f ` set_of X)
  using compact_continuous_image[of "set_of X" f, simplified assms] compact_interval
  unfolding set_of_eq
  by simp

lemma sup_inf_dist_bounded:
  fixes X::real set
  shows bdd_below X  bdd_above X ==>   x  X.  x'  X. dist x  x'  Sup X - Inf X
  using cInf_lower[of _ X] cSup_upper[of _ X]
  apply(auto simp add:dist_real_def)[1]
  by (smt (z3))

lemma set_of_nonempty[simp]: set_of X  {}
  by(simp add:set_of_eq)

lemma lower_in_interval[simp]:lower X i X
  by(simp add: in_intervalI)

lemma upper_in_interval[simp]:upper X i X
  by(simp add: in_intervalI)

lemma bdd_below_set_of: bdd_below (set_of X)
  by (metis atLeastAtMost_iff bdd_below.unfold set_of_eq)

lemma bdd_above_set_of: bdd_above (set_of X)
  by (metis atLeastAtMost_iff bdd_above.unfold set_of_eq)

lemma closed_set_of: closed (set_of (X::real interval))
  by (metis  closed_real_atLeastAtMost set_of_eq)

lemma set_f_nonempty: f ` set_of X  {}
  apply(simp add:image_def set_of_eq)
  by fastforce

lemma interval_linorder_case_split[case_names LeftOf Including  RightOf]:
assumes (upper x < c  P (x::('a::linorder interval)))
   ( c i x  P x)
   ( c < lower x  P x)
 showsP x
proof(insert assms, cases " upper x < c")
  case True
  then show ?thesis
    using assms(1) by blast
next
  case False
  then show ?thesis
  proof(cases "c i x")
    case True
    then show ?thesis
      using assms(2) by blast
  next
    case False
    then show ?thesis
      by (meson assms(1) assms(3) in_intervalI not_le_imp_less)
  qed
qed

lemma foldl_conj_True:
foldl (∧) x xs = list_all (λ e. e = True)  (x#xs)
  by(induction xs rule:rev_induct, auto)

lemma foldl_conj_set_True:
foldl (∧) x xs = ( e  set (x#xs) . e = True)
  by(induction xs rule:rev_induct, auto)


section "Interval Bounds and Set Conversion "

lemma sup_set_of:
  fixes X :: "'a::{conditionally_complete_lattice} interval"
  shows "Sup (set_of X) = upper X"
  unfolding set_of_def upper_def
  using cSup_atLeastAtMost lower_le_upper[of X]
  by (simp add: upper_def lower_def)

lemma inf_set_of:
  fixes X :: "'a::{conditionally_complete_lattice} interval"
  shows "Inf (set_of X) = lower X"
  unfolding set_of_def lower_def
  using cInf_atLeastAtMost lower_le_upper[of X]
  by (simp add: upper_def lower_def)

lemma inf_le_sup_set_of:
  fixes X :: "'a::{conditionally_complete_lattice} interval"
  shows"Inf (set_of X)  Sup (set_of X)"
  using sup_set_of inf_set_of lower_le_upper
  by metis

lemma in_bounds: x i X  lower X  x  x  upper X
  by (simp add: set_of_eq)

lemma lower_bounds[simp]:
  assumes L  U
  shows lower (Interval(L,U)) = L
  using assms
  apply (simp add: lower.rep_eq)
  by (simp add: bounds_of_interval_eq_lower_upper lower_Interval upper_Interval)

lemma upper_bounds [simp]:
  assumes L  U
  shows upper (Interval(L,U)) = U
  using assms
  apply (simp add: upper.rep_eq)
  by (simp add: bounds_of_interval_eq_lower_upper lower_Interval upper_Interval)

lemma lower_point_interval[simp]: lower (Interval (x,x)) = x
  by (simp)
lemma upper_point_interval[simp]: upper (Interval (x,x)) = x
  by (simp)

lemma map2_nth:
  assumes  length xs = length ys
      and      n < length xs
    shows (map2 f xs ys)!n = f (xs!n) (ys!n)
  using assms  by simp

lemma map_set: a  set (map f X)   ( x  set X . f x = a)
  by auto
lemma map_pair_set_left: (a,b)  set (zip (map f X) (map f Y))  ( x  set X . f x = a)
  by (meson map_set set_zip_leftD)
lemma map_pair_set_right: (a,b)  set (zip (map f X) (map f Y))   ( y  set Y . f y = b)
  by (meson map_set set_zip_rightD)
lemma map_pair_set: (a,b)  set (zip (map f X) (map f Y))   ( x  set X . f x = a)   ( y  set Y . f y = b)
  by (meson map_pair_set_left set_zip_rightD zip_same)

lemma map_pair_f_all:
  assumes length X = length Y
  shows ((x, y)set (zip (map f X) (map f Y)). x  y) = ((x, y)set (zip X Y ). f x  f y)
  by(insert assms(1), induction X Y rule:list_induct2, auto)

definition   map_interval_swap :: ('a::linorder × 'a) list  'a interval list where
map_interval_swap = map (λ (x,y). Interval (if x  y then (x,y) else (y,x)))

definition  mk_interval :: ('a::linorder × 'a)  'a interval where
mk_interval = (λ (x,y). Interval (if x  y then (x,y) else (y,x)))

definition  mk_interval' :: ('a::linorder × 'a)  'a interval where
mk_interval' = (λ (x,y). (if x  y then Interval(x,y) else Interval(y,x)))

lemma map_interval_swap_code[code]:
  map_interval_swap = map (λ (x,y). the (if x  y then Interval' x y else Interval' y x))
  unfolding map_interval_swap_def by(rule ext, rule arg_cong, auto simp add: Interval'.abs_eq)

lemma  mk_interval_code[code]:
  mk_interval = (λ (x,y). the (if x  y then Interval' x y else Interval' y x))
  unfolding mk_interval_def by(rule ext, rule arg_cong, auto simp add: Interval'.abs_eq)

lemma  mk_interval':
  mk_interval = (λ (x,y). (if x  y then Interval(x, y) else Interval(y, x)))
  unfolding mk_interval_def by(rule ext, rule arg_cong, auto)

lemma mk_interval_lower[simp]: lower (mk_interval (x,y)) = (if x  y then x else y)
  by(simp add: lower_def Interval_inverse mk_interval_def)

lemma mk_interval_upper[simp]: upper (mk_interval (x,y)) = (if x  y then y else x)
  by(simp add: upper_def Interval_inverse mk_interval_def)

section‹Linear Order on List of Intervals›

definition
  le_interval_list :: ('a::linorder) interval list  'a interval list  bool ("(_/ I _)" [51, 51] 50)
  where
    le_interval_list Xs Ys  (length Xs = length Ys)  (foldl (∧) True (map2 (≤) Xs Ys))

lemma le_interval_single: (x  y) = ([x] I [y])
  unfolding le_interval_list_def by simp

lemma le_intervall_empty[simp]: [] I []
  unfolding le_interval_list_def by simp

lemma le_interval_list_rev: (is I js) = (rev is I rev js)
  unfolding le_interval_list_def
  by(safe, simp_all add: foldl_conj_set_True zip_rev)

lemma le_interval_list_imp_length:
  assumes Xs I Ys shows length Xs = length Ys
using assms unfolding le_interval_list_def
  by simp


lemma lsplit_left: assumes length (xs) = length (ys)
  and (n<length (x # xs). (x # xs) ! n  (y # ys) ! n) shows (( n < length xs. xs!n  ys!n)  x  y)
  using assms  by auto

lemma lsplit_right: assumes length (xs) = length (ys)
  and  (( n < length xs. xs!n  ys!n)  x  y)
shows (n<length (x # xs)  (x # xs) ! n  (y # ys) ! n)
proof(cases "n=0")
  case True
  then show ?thesis using assms by(simp)
next
  case False
  then show ?thesis using assms by simp
qed

lemma lsplit: assumes length (xs) = length (ys)
  shows (n<length (x # xs). (x # xs) ! n  (y # ys) ! n) =
      (( n < length xs. xs!n  ys!n)  x  y)
  using assms lsplit_left lsplit_right by metis

lemma le_interval_list_all':
  assumes length Xs = length Ys and Xs I Ys shows  n < length Xs. Xs!n  Ys!n
proof(insert assms, induction rule:list_induct2)
  case Nil
  then show ?case by simp
next
  case (Cons x xs y ys)
  then show ?case
    using lsplit[of xs ys x y] le_interval_list_def
    unfolding le_interval_list_def
    by(auto simp add: foldl_conj_True)
qed

lemma le_interval_list_all2:
  assumes length Xs = length Ys
    and  n<length Xs . (Xs!n   Ys!n)
  shows  Xs I Ys
proof(insert assms, induction rule:list_induct2)
  case Nil
  then show ?case by simp
next
  case (Cons x xs y ys)
  then show ?case
    using lsplit[of xs ys x y] le_interval_list_def
    unfolding le_interval_list_def
    by(auto simp add: foldl_conj_True)
qed

lemma le_interval_list_all:
  assumes Xs I Ys shows  n < length Xs. Xs!n  Ys!n
  using assms le_interval_list_all' le_interval_list_imp_length
  by auto

lemma le_interval_list_imp:
  assumes Xs I Ys shows n < length Xs   Xs!n  Ys!n
  using assms le_interval_list_all' le_interval_list_imp_length
  by auto

lemma  interval_set_leq_eq: (X  Y) = (set_of X  set_of Y)
  for X :: 'a::linordered_ring interval
  by (simp add: less_eq_interval_def set_of_eq)

lemma times_interval_right:
  fixes X Y C ::'a::linordered_ring interval
  assumes X  Y
  shows C * X  C * Y
  using assms[simplified interval_set_leq_eq]
  apply(subst interval_set_leq_eq)
  by (simp add: set_of_mul_inc_right)


lemma times_interval_left:
  fixes X Y C ::'a::{real_normed_algebra,linordered_ring,linear_continuum_topology} interval
  assumes X  Y
  shows X * C  Y * C
  using assms[simplified interval_set_leq_eq]
  apply(subst interval_set_leq_eq)
  by (simp add: set_of_mul_inc_left)

section‹Support for Lists of Intervals ›

abbreviation in_interval_list::('a::preorder) list  'a interval list  bool ("(_/ I _)" [51, 51] 50)
  where in_interval_list xs Xs  foldl (∧) True (map2 (in_interval) xs Xs)

lemma interval_of_in_interval_list[simp]: xs I map interval_of xs
proof(induction xs)
  case Nil
  then show ?case by simp
next
  case (Cons a xs)
  then show ?case
    by (simp add: in_interval_to_interval)
qed

lemma interval_of_in_eq:  interval_of x  X = (x i X)
  by (simp add: less_eq_interval_def set_of_eq)

lemma interval_of_list_in:
  assumes length inputs = length Inputs
shows (map interval_of inputs I Inputs) = (inputs I Inputs)
unfolding le_interval_list_def
proof(insert assms, induction inputs Inputs rule:list_induct2)
  case Nil
  then show ?case by simp
next
  case (Cons x xs y ys)
  then show ?case
  proof(cases "interval_of x  y")
    case True
    then show ?thesis
      using Cons by(simp add: interval_of_in_eq)
  next
    case False
    then show ?thesis using foldl_conj_True interval_of_in_eq by auto
  qed
qed

section ‹Interval Width and Arithmetic Operations›

lemma interval_width_addition:
  fixes A:: "'a::{linordered_ring} interval"
  shows width (A + B) = width A + width B
  by (simp add: width_def)

lemma interval_width_times:
  fixes a :: "'a::{linordered_ring}" and A :: "'a interval"
  shows "width (interval_of a * A) = ¦a¦ * width A"
proof -
  have "width (interval_of a * A) = upper (interval_of a * A) - lower (interval_of a * A)"
    by (simp add: width_def)
  also have "... = (if a > 0 then a * upper A - a * lower A else a * lower A - a * upper A)"
    proof -
      have "upper (interval_of a * A) = max (a * lower A) (a * upper A)"
        by (simp add: upper_times)
      moreover have "lower (interval_of a * A) = min (a * lower A) (a * upper A)"
        by (simp add: lower_times)
      ultimately show ?thesis
        by (simp add: mult_left_mono mult_left_mono_neg)
    qed
  also have "... = ¦a¦ * (upper A - lower A)"
    by (simp add: right_diff_distrib)
  finally show ?thesis
    by (simp add: width_def)
qed

lemma interval_sup_width:
  fixes X Y :: "'a::{linordered_ring, lattice} interval"
  shows "width (sup X Y) = max (upper X) (upper Y) - min (lower X) (lower Y)"
proof -
  have a0: "i ia. lower (sup i ia) = min (lower i::real) (lower ia)"
    by (simp add: inf_min)
  have a1: "i ia. upper (sup i ia) = max (upper i::real) (upper ia)"
    by (simp add: sup_max)
  show ?thesis
    using a0 a1 by (simp add: inf_min sup_max width_def)
qed

lemma width_expanded: interval_of (width Y) = Interval(upper Y - lower Y, upper Y - lower Y)
  unfolding width_def interval_of_def by simp

lemma interval_width_positive:
  fixes X :: 'a::{linordered_ring} interval
  shows 0  width X
  using lower_le_upper
  by (simp add: width_def)

section ‹Interval Multiplication›

lemma interval_interval_times:
X * Y = Interval(Min {(lower X * lower Y), (lower X * upper Y), (upper X * lower Y), (upper X * upper Y)},
                  Max {(lower X * lower Y), (lower X * upper Y), (upper X * lower Y), (upper X * upper Y)})
  by (simp add: times_interval_def bounds_of_interval_eq_lower_upper Let_def min_def)

lemma interval_times_scalar: interval_of a * A = mk_interval(a * lower A, a * upper A)
 proof -
    have "upper (interval_of a * A) = max (a * lower A) (a * upper A)"
      by (simp add: upper_times )
    moreover have "lower (interval_of a * A) = min (a * lower A) (a * upper A)"
      by (simp add: lower_times)
    ultimately show ?thesis
      by (simp add: interval_eq_iff)
  qed

lemma interval_times_scalar_pos_l:
  fixes a :: "'a::{ordered_semiring,times,linorder}"
  assumes 0  a
  shows interval_of a * A = Interval(a * lower A, a * upper A)
 proof -
   have "upper (interval_of a * A) = a * upper A"
     using assms
     by (simp add: upper_times mult_left_mono)
   moreover have "lower (interval_of a * A) = a * lower A"
     using assms
     by (simp add: lower_times mult_left_mono)
   ultimately show ?thesis
     using assms
     by (metis bounds_of_interval_eq_lower_upper bounds_of_interval_inverse lower_le_upper)
 qed

lemma interval_times_scalar_pos_r:
  fixes a :: "'a::{linordered_idom}"
  assumes 0  a
  shows A * interval_of a = Interval(a * lower A, a * upper A)
 proof -
   have "upper (A * interval_of a) = a * upper A"
     using assms
     by (simp add: mult.commute mult_left_mono upper_times)
   moreover have "lower (A * interval_of a ) = a * lower A"
     using assms
     by (simp add: lower_times mult_right_mono)
   ultimately show ?thesis
     using assms
     by (metis bounds_of_interval_eq_lower_upper bounds_of_interval_inverse lower_le_upper)
 qed


section ‹Distance-based Properties of Intervals›

text ‹Given two real intervals @{term X} and @{term Y}, and two real numbers  @{term a}
and @{term b}, the width of the sum of the scaled intervals is equivalent to the width of the
two individual intervals.›

lemma width_of_scaled_interval_sum:
  fixes  X :: "'a::{linordered_ring} interval"
  shows width (interval_of a * X + interval_of b * Y) = ¦a¦ * width X + ¦b¦ * width Y
  using interval_width_addition interval_width_times by metis

lemma width_of_product_interval_bound_real:
  fixes X :: "real interval"
  shows interval_of (width (X * Y))  abs_interval(X) * interval_of (width Y) + abs_interval(Y) * interval_of (width X)
proof -
  have a0: "lower X  upper X" lower Y  upper Y
    using lower_le_upper by simp_all
  have a1: width Y  0 and a1': width X  0
    using a0 interval_width_positive by simp_all
  have a2: abs_interval(X) * interval_of (width Y) = Interval (width Y * lower (abs_interval X), width Y * upper (abs_interval X))
    using interval_times_scalar_pos_r a0 a1 by blast
  have  a3: abs_interval(Y) * interval_of (width X) = Interval (width X * lower (abs_interval Y), width X * upper (abs_interval Y))
    using interval_times_scalar_pos_r a0 interval_width_positive by blast
  have a4: abs_interval(X) * interval_of (width Y) + abs_interval(Y) * interval_of (width X) = Interval (width Y * lower (abs_interval X) + width X * lower (abs_interval Y), width Y * upper (abs_interval X) + width X * upper (abs_interval Y))
    using a0 a2 a3 unfolding plus_interval_def
    apply (simp add:bounds_of_interval_eq_lower_upper mult_left_mono interval_width_positive)
    by (auto simp add:bounds_of_interval_eq_lower_upper mult_left_mono interval_width_positive)
  have a5: width(X * Y) = Max {(lower X * lower Y), (lower X * upper Y), (upper X * lower Y), (upper X * upper Y)} - Min {(lower X * lower Y), (lower X * upper Y), (upper X * lower Y), (upper X * upper Y)}
    using lower_times upper_times unfolding  width_def by metis
  have a6: Max {(lower X * lower Y), (lower X * upper Y), (upper X * lower Y), (upper X * upper Y)} - Min {(lower X * lower Y), (lower X * upper Y), (upper X * lower Y), (upper X * upper Y)}  width Y * upper (abs_interval X) + width X * upper (abs_interval Y)
    using a0 a1 a1' unfolding width_def
    by (simp, smt (verit) a0  mult.commute mult_less_cancel_left_pos mult_minus_left right_diff_distrib)  (* takes long *)
    have a7: width Y * lower (abs_interval X) + width X * lower (abs_interval Y)  Max {(lower X * lower Y), (lower X * upper Y), (upper X * lower Y), (upper X * upper Y)} - Min {(lower X * lower Y), (lower X * upper Y), (upper X * lower Y), (upper X * upper Y)}
    using a0 a1 a1' unfolding width_def
    by (simp, smt (verit) a0 mult.commute mult_less_cancel_left_pos mult_minus_left right_diff_distrib) (* takes long *)
  have interval_of (Max {(lower X * lower Y), (lower X * upper Y), (upper X * lower Y), (upper X * upper Y)} -
                     Min {(lower X * lower Y), (lower X * upper Y), (upper X * lower Y), (upper X * upper Y)}) 
        Interval (width Y * lower (abs_interval X) + width X * lower (abs_interval Y),
                  width Y * upper (abs_interval X) + width X * upper (abs_interval Y))
    using a6 a7 unfolding less_eq_interval_def
    by (smt (verit, ccfv_threshold) lower_bounds lower_interval_of upper_bounds upper_interval_of)
  then show ?thesis using a4 a5 a6 by presburger
qed


lemma width_of_product_interval_bound_int:
  fixes X :: "int interval"
  shows interval_of (width (X * Y))  abs_interval(X) * interval_of (width Y) + abs_interval(Y) * interval_of (width X)
proof -
  have a0: "lower X  upper X" lower Y  upper Y
    using lower_le_upper by simp_all
  have a1: width Y  0 and width X  0
    using a0 interval_width_positive by simp_all
   have a2: abs_interval(X) * interval_of (width Y) = Interval (width Y * lower (abs_interval X), width Y * upper (abs_interval X))
    using interval_times_scalar_pos_r a0 a1 by blast
  have  a3: abs_interval(Y) * interval_of (width X) = Interval (width X * lower (abs_interval Y), width X * upper (abs_interval Y))
    using interval_times_scalar_pos_r a0 interval_width_positive by blast
  have a4: abs_interval(X) * interval_of (width Y) + abs_interval(Y) * interval_of (width X) = Interval (width Y * lower (abs_interval X) + width X * lower (abs_interval Y), width Y * upper (abs_interval X) + width X * upper (abs_interval Y))
    using a0 a2 a3 unfolding plus_interval_def
    by (auto simp add:bounds_of_interval_eq_lower_upper mult_left_mono interval_width_positive) (* takes long *)
  have a5: width(X * Y) = Max {(lower X * lower Y), (lower X * upper Y), (upper X * lower Y), (upper X * upper Y)} - Min {(lower X * lower Y), (lower X * upper Y), (upper X * lower Y), (upper X * upper Y)}
    using lower_times upper_times unfolding  width_def by metis
  have a6: Max {(lower X * lower Y), (lower X * upper Y), (upper X * lower Y), (upper X * upper Y)} - Min {(lower X * lower Y), (lower X * upper Y), (upper X * lower Y), (upper X * upper Y)}  width Y * upper (abs_interval X) + width X * upper (abs_interval Y)
    using a0 unfolding width_def
    by (simp, smt (verit) a0  mult.commute mult_less_cancel_left_pos mult_minus_left right_diff_distrib)  (* takes long *)
  have a7: width Y * lower (abs_interval X) + width X * lower (abs_interval Y)  Max {(lower X * lower Y), (lower X * upper Y), (upper X * lower Y), (upper X * upper Y)} - Min {(lower X * lower Y), (lower X * upper Y), (upper X * lower Y), (upper X * upper Y)}
    using a0 unfolding width_def
    apply(auto)[1]
    by (smt (verit) a0(1) a0(2) left_diff_distrib mult_less_cancel_left_pos int_distrib(3)
        mult_less_cancel_left mult_minus_right mult.commute mult_le_0_iff right_diff_distrib'
        mult_left_mono mult_le_cancel_right right_diff_distrib zero_less_mult_iff )+
  have interval_of (Max {(lower X * lower Y), (lower X * upper Y), (upper X * lower Y), (upper X * upper Y)} -
                     Min {(lower X * lower Y), (lower X * upper Y), (upper X * lower Y), (upper X * upper Y)}) 
        Interval (width Y * lower (abs_interval X) + width X * lower (abs_interval Y),
                  width Y * upper (abs_interval X) + width X * upper (abs_interval Y))
    using a6 a7 unfolding less_eq_interval_def
    by (smt (verit, ccfv_threshold) lower_bounds lower_interval_of upper_bounds upper_interval_of)
  then show ?thesis using a4 a5 a6 by presburger
qed

end