From 5bf5642ecebcf1bb870520e6f5441c64cc313df9 Mon Sep 17 00:00:00 2001 From: lecopivo Date: Thu, 25 Jul 2024 17:16:14 -0400 Subject: [PATCH] cleaning up param deriv examples --- .../FunctionTransformations/Preimage.lean | 18 +- SciLean/Core/Functions/Gaussian.lean | 62 +++++ .../HasParamDerivWithJumps/Common.lean | 41 +++- .../HasParamDerivWithJumps/Functions.lean | 6 +- .../HasParamFDerivWithJumps.lean | 66 ++--- .../HasParamFwdFDerivWithJumps.lean | 41 ++-- .../HasParamRevFDerivWithJumps.lean | 60 ++--- SciLean/Doodle.lean | 228 +++--------------- test/param_deriv/affine_discont_2d.lean | 4 +- test/param_deriv/circle_discont_2d.lean | 6 +- test/param_deriv/risk_minimization.lean | 33 ++- test/param_deriv/simple_discont_1d.lean | 2 +- test/param_deriv/simple_discont_2d.lean | 4 +- 13 files changed, 255 insertions(+), 316 deletions(-) diff --git a/SciLean/Core/FunctionTransformations/Preimage.lean b/SciLean/Core/FunctionTransformations/Preimage.lean index 02fc6a5f..d08da4f4 100644 --- a/SciLean/Core/FunctionTransformations/Preimage.lean +++ b/SciLean/Core/FunctionTransformations/Preimage.lean @@ -38,53 +38,53 @@ theorem preimage_comp' (f : β → γ) (g : α → β) : theorem Prod.mk.arg_fstsnd.preimage_rule_prod (f : α → β) (g : α → γ) (B : Set β) (C : Set γ) : preimage (fun x => (f x, g x)) (B.prod C) = - f ⁻¹' B ∩ g ⁻¹' C := sorry_proof + f ⁻¹' B ∩ g ⁻¹' C := by ext; simp[Set.prod] @[fun_trans] theorem Prod.mk.arg_fst.preimage_rule_prod (f : α → β) (c : γ) : preimage (fun x => (f x, c)) = - fun s => f ⁻¹' (s.fst c) := sorry_proof + fun s => f ⁻¹' (s.fst c) := by ext; simp @[fun_trans] theorem Prod.mk.arg_snd.preimage_rule_prod (b : β) (g : α → γ) : preimage (fun x => (b, g x)) = - fun s => g ⁻¹' (s.snd b) := sorry_proof + fun s => g ⁻¹' (s.snd b) := by ext; simp open SciLean variable {R} [RealScalar R] -- probably generalize following to LinearlyOrderedAddCommGroup @[fun_trans] -theorem HAdd.hAdd.arg_a0.preimage_rule_Ioo (x' a b : R) : +theorem HAdd.hAdd.arg_a0.preimage_rule_Ioo (x' a b : R) : preimage (fun x : R => x + x') (Ioo a b) = - Ioo (a - x') (b - x') := by ext; simp; sorry_proof + Ioo (a - x') (b - x') := by ext; simp; sorry_proof -- over ℝ `simp` finishes the proof @[fun_trans] theorem HAdd.hAdd.arg_a1.preimage_rule_Ioo (x' a b : R) : preimage (fun x : R => x' + x) (Ioo a b) = - Ioo (a - x') (b - x') := by ext; simp; sorry_proof + Ioo (a - x') (b - x') := by ext; simp; sorry_proof -- over ℝ `simp` finishes the proof @[fun_trans] theorem HSub.hSub.arg_a0.preimage_rule_Ioo (x' a b : R) : preimage (fun x : R => x - x') (Ioo a b) = - Ioo (a + x') (b + x') := by ext; simp; sorry_proof + Ioo (a + x') (b + x') := by ext; simp; sorry_proof -- over ℝ `simp` finishes the proof @[fun_trans] theorem HSub.hSub.arg_a1.preimage_rule_Ioo (x' a b : R) : preimage (fun x : R => x' - x) (Ioo a b) = - Ioo (x' - b) (x' - a) := by ext; simp; sorry_proof + Ioo (x' - b) (x' - a) := by ext; simp; sorry_proof -- over ℝ `simp` finishes the proof @[fun_trans] theorem Neg.neg.arg_a1.preimage_rule_Ioo (a b : R) : preimage (fun x : R => - x) (Ioo a b) = - Ioo (-b) (-a) := by ext; simp; sorry_proof + Ioo (-b) (-a) := by ext; simp; sorry_proof -- over ℝ `simp` finishes the proof diff --git a/SciLean/Core/Functions/Gaussian.lean b/SciLean/Core/Functions/Gaussian.lean index 9c8b54d6..ec2340ad 100644 --- a/SciLean/Core/Functions/Gaussian.lean +++ b/SciLean/Core/Functions/Gaussian.lean @@ -1,5 +1,6 @@ import SciLean.Core.FunctionTransformations import SciLean.Core.Functions.Exp +import SciLean.Core.Functions.Inner import SciLean.Tactic.Autodiff -- import SciLean.Core.Meta.GenerateRevDeriv @@ -69,6 +70,67 @@ theorem gaussian.arg_μx.Differentiable_rule Differentiable R (fun w => gaussian (μ w) σ (x w)) := by intro w; fun_prop +set_option linter.unusedVariables false in +@[fun_trans] +theorem gaussian.arg_μx.fderiv_rule + {W} [NormedAddCommGroup W] [NormedSpace R W] + {U} [NormedAddCommGroup U] [AdjointSpace R U] + (μ : W → U) (σ : R) (x : W → U) + (hμ : Differentiable R μ) (hx : Differentiable R x) : + fderiv R (fun w => gaussian (μ w) σ (x w)) + = + fun w => fun dw =>L[R] + let μ' := μ w + let dμ := fderiv R μ w dw + let x' := x w + let dx := fderiv R x w dw + let dx' := - σ^(-2:ℤ) * ⟪dx-dμ, x'-μ'⟫ + dx' * gaussian μ' σ x' := sorry_proof + + +set_option linter.unusedVariables false in +@[fun_trans] +theorem gaussian.arg_μx.fwdFDeriv_rule + {W} [NormedAddCommGroup W] [NormedSpace R W] + {U} [NormedAddCommGroup U] [AdjointSpace R U] + (μ : W → U) (σ : R) (x : W → U) + (hμ : Differentiable R μ) (hx : Differentiable R x) : + fwdFDeriv R (fun w => gaussian (μ w) σ (x w)) + = + fun w dw => + let μdμ := fwdFDeriv R μ w dw + let xdx := fwdFDeriv R x w dw + let dx' := - σ^(-2:ℤ) * ⟪xdx.2-μdμ.2, xdx.1-μdμ.1⟫ + let s := gaussian μdμ.1 σ xdx.1 + (s, dx' * s) := by unfold fwdFDeriv; autodiff + + +set_option linter.unusedVariables false in +@[fun_trans] +theorem gaussian.arg_μx.revFDeriv_rule + {W} [NormedAddCommGroup W] [AdjointSpace R W] [CompleteSpace W] + {U} [NormedAddCommGroup U] [AdjointSpace R U] [CompleteSpace U] + (μ : W → U) (σ : R) (x : W → U) + (hμ : Differentiable R μ) (hx : Differentiable R x) : + revFDeriv R (fun w => gaussian (μ w) σ (x w)) + = + fun w => + let μdμ := revFDeriv R μ w + let xdx := revFDeriv R x w + let s := gaussian μdμ.1 σ xdx.1 + (s, + fun dr => + let dx' := (dr * s * σ^(-2:ℤ)) • (xdx.1-μdμ.1) + let dw₁ := xdx.2 dx' + let dw₂ := μdμ.2 dx' + dw₂ - dw₁) := by + unfold revFDeriv + funext w; simp + autodiff; autodiff + funext dw; simp[← smul_assoc,mul_comm] + sorry_proof + + @[fun_prop] theorem gaussian.arg_μx.CDifferentiableAt_rule (w : W) (μ : W → U) (σ : R) (x : W → U) diff --git a/SciLean/Core/Transformations/HasParamDerivWithJumps/Common.lean b/SciLean/Core/Transformations/HasParamDerivWithJumps/Common.lean index ac089eb3..f95454a3 100644 --- a/SciLean/Core/Transformations/HasParamDerivWithJumps/Common.lean +++ b/SciLean/Core/Transformations/HasParamDerivWithJumps/Common.lean @@ -130,7 +130,7 @@ set_default_scalar R -- frontierSpeed is not well defined @[simp,ftrans_simp] theorem frontierSpeed_setOf_le (φ ψ : W → X → R) : - frontierSpeed' R (fun w => {x | φ w x ≤ ψ w x}) + frontierSpeed R (fun w => {x | φ w x ≤ ψ w x}) = fun w dw x => let ζ := (fun w x => φ w x - ψ w x) @@ -139,7 +139,7 @@ theorem frontierSpeed_setOf_le (φ ψ : W → X → R) : @[simp,ftrans_simp] theorem frontierSpeed_setOf_lt (φ ψ : W → X → R) : - frontierSpeed' R (fun w => {x | φ w x < ψ w x}) + frontierSpeed R (fun w => {x | φ w x < ψ w x}) = fun w dw x => let ζ := (fun w x => φ w x - ψ w x) @@ -147,6 +147,43 @@ theorem frontierSpeed_setOf_lt (φ ψ : W → X → R) : sorry_proof +-- not sure what to do when `(a w) > (b w)`. In that case is not really well defined `frontierSpeed` +set_option linter.unusedVariables false in +open Set in +@[simp,ftrans_simp] +theorem frontierSpeed_Icc (a b : W → R) (ha : Differentiable R a) (hb : Differentiable R b) : + frontierSpeed R (fun w => Icc (a w) (b w)) + = + fun w dw x => + let ada := fwdFDeriv R a w dw + let bdb := fwdFDeriv R b w dw + let m := (ada.1 + bdb.1)/2 + if x < m then + - ada.2 + else + bdb.2 := by + sorry_proof + + +set_option linter.unusedVariables false in +open Set in +@[simp,ftrans_simp] +theorem frontierGrad_Icc + {W} [NormedAddCommGroup W] [AdjointSpace R W] [CompleteSpace W] + (a b : W → R) (ha : Differentiable R a) (hb : Differentiable R b) : + frontierGrad R (fun w => Icc (a w) (b w)) + = + fun w x => + let ada := revFDeriv R a w + let bdb := revFDeriv R b w + let m := (ada.1 + bdb.1)/2 + if x < m then + - ada.2 1 + else + bdb.2 1 := by + sorry_proof + + end ---------------------------------------------------------------------------------------------------- diff --git a/SciLean/Core/Transformations/HasParamDerivWithJumps/Functions.lean b/SciLean/Core/Transformations/HasParamDerivWithJumps/Functions.lean index 58c3fad6..50ee606b 100644 --- a/SciLean/Core/Transformations/HasParamDerivWithJumps/Functions.lean +++ b/SciLean/Core/Transformations/HasParamDerivWithJumps/Functions.lean @@ -23,7 +23,7 @@ variable open Scalar in @[gtrans] def Norm2.norm2.HasParamFDerivWithJumpsAt_rule := - (HasParamFDerivWithJumpsAt.comp1_smooth_jumps_rule (R:=R) (W:=W) (X:=X) (Y:=U) (Z:=R) + (HasParamFDerivWithJumpsAt.comp1_differentiable_jumps_rule (R:=R) (W:=W) (X:=X) (Y:=U) (Z:=R) (fun _ y => ‖y‖₂²[R]) (by simp; fun_prop)) rewrite_type_by (repeat ext); autodiff @@ -31,14 +31,14 @@ def Norm2.norm2.HasParamFDerivWithJumpsAt_rule := open Scalar in @[gtrans] def Norm2.norm2.HasParamFwdFDerivWithJumpsAt_rule := - (HasParamFwdDerivWithJumps.comp1_smooth_jumps_rule (R:=R) (W:=W) (X:=X) (Y:=U) (Z:=R) + (HasParamFwdDerivWithJumps.comp1_differentiable_jumps_rule (R:=R) (W:=W) (X:=X) (Y:=U) (Z:=R) (fun _ y => ‖y‖₂²[R]) (by simp; fun_prop)) rewrite_type_by (repeat ext); autodiff open Scalar in @[gtrans] def Norm2.norm2.HasParamRevFDerivWithJumpsAt_rule := - (HasParamRevFDerivWithJumpsAt.comp1_smooth_jumps_rule (R:=R) (W:=U) (X:=X) (Y:=V) (Z:=R) + (HasParamRevFDerivWithJumpsAt.comp1_differentiable_jumps_rule (R:=R) (W:=U) (X:=X) (Y:=V) (Z:=R) (fun _ y => ‖y‖₂²[R]) (by simp; fun_prop)) rewrite_type_by (repeat ext); autodiff diff --git a/SciLean/Core/Transformations/HasParamDerivWithJumps/HasParamFDerivWithJumps.lean b/SciLean/Core/Transformations/HasParamDerivWithJumps/HasParamFDerivWithJumps.lean index 2ecba7e5..7eb0adf7 100644 --- a/SciLean/Core/Transformations/HasParamDerivWithJumps/HasParamFDerivWithJumps.lean +++ b/SciLean/Core/Transformations/HasParamDerivWithJumps/HasParamFDerivWithJumps.lean @@ -31,7 +31,7 @@ set_default_scalar R variable (R) open Classical in noncomputable -def frontierSpeed' (A : W → Set X) (w dw : W) (x : X) : R := +def frontierSpeed (A : W → Set X) (w dw : W) (x : X) : R := match Classical.dec (∃ (φ : W → X → R), (∀ w, closure (A w) = {x | φ w x ≤ 0})) with | .isTrue h => let φ := Classical.choose h @@ -79,11 +79,11 @@ structure HasParamFDerivWithJumpsAtImpl (f : W → X → Y) (w : W) ∀ p n : J, match ι p n with | none => True | some i => ∀ x ∈ jump i, - frontierSpeed' R (Ω n) w dw x = jumpSpeed i dw x + frontierSpeed R (Ω n) w dw x = jumpSpeed i dw x variable {R} variable (R W X Y) -structure DiscontinuityData where +structure DiscontinuityData (R W X Y : Type) where vals : X → Y×Y speed : W → X → R discontinuity : Set X @@ -141,28 +141,28 @@ theorem fderiv_under_integral = let interior := ∫ x, f' dw x ∂μ let density := fun x => Scalar.ofENNReal (R:=R) (μ.rnDeriv volume x) - let shocks := disc.foldl (init:=0) - fun sum ⟨df,s,S⟩ => sum + ∫ x in S, - let vals := df x - (s dw x * density x) • (vals.1 - vals.2) ∂μH[finrank R X - (1:ℕ)] + let shocks := disc.foldl (init:=0) fun sum ⟨df,s,Γ⟩ => sum + + ∫ x in Γ, + let vals := df x + (s dw x * density x) • (vals.1 - vals.2) ∂μH[finrank R X - (1:ℕ)] interior + shocks := sorry_proof -- @[fun_trans] theorem fderiv_under_integral_over_set {X} [NormedAddCommGroup X] [AdjointSpace R X] [NormedSpace ℝ X] [CompleteSpace X] [MeasureSpace X] [BorelSpace X] - (f : W → X → Y) (w dw : W) (μ : Measure X) (A : Set X) + (f : W → X → Y) (w dw : W) (μ : Measure X) (Ω : Set X) {f' disc} (hf : HasParamFDerivWithJumpsAt R f w f' disc) - (hA : AlmostDisjoint (frontier A) disc.getDiscontinuity μH[finrank ℝ X - (1:ℕ)]) + (hA : AlmostDisjoint (frontier Ω) disc.getDiscontinuity μH[finrank ℝ X - (1:ℕ)]) /- todo: add some integrability conditions -/ : - (fderiv R (fun w' => ∫ x in A, f w' x ∂μ) w dw) + (fderiv R (fun w' => ∫ x in Ω, f w' x ∂μ) w dw) = - let interior := ∫ x in A, f' dw x ∂μ + let interior := ∫ x in Ω, f' dw x ∂μ let density := fun x => Scalar.ofENNReal (R:=R) (μ.rnDeriv volume x) let shocks := disc.foldl (init:=0) - fun sum ⟨df,s,S⟩ => sum + - ∫ x in S ∩ A, + fun sum ⟨df,s,Γ⟩ => sum + + ∫ x in Γ ∩ Ω, let vals := df x (s dw x * density x) • (vals.1 - vals.2) ∂μH[finrank R X - (1:ℕ)] interior + shocks := sorry_proof @@ -177,7 +177,7 @@ namespace HasParamFDerivWithJumpsAt @[gtrans high] -theorem smooth_rule +theorem differentiable_at_rule (w : W) (f : W → X → Y) (hf : ∀ x, DifferentiableAt R (f · x) w) : HasParamFDerivWithJumpsAt R f w @@ -186,7 +186,7 @@ theorem smooth_rule sorry_proof -theorem comp_smooth_jumps_rule_at +theorem comp_differentiable_jumps_rule_at (f : W → Y → Z) (g : W → X → Y) (w : W) {g' disc} (hf : ∀ x, DifferentiableAt R (fun (w,y) => f w y) (w,g w x)) @@ -207,7 +207,7 @@ theorem comp_smooth_jumps_rule_at -theorem comp_smooth_jumps_rule +theorem comp_differentiable_jumps_rule (f : W → Y → Z) (g : W → X → Y) (w : W) {g' disc} (hf : Differentiable R (fun (w,y) => f w y)) @@ -228,7 +228,7 @@ theorem comp_smooth_jumps_rule -theorem comp1_smooth_jumps_rule +theorem comp1_differentiable_jumps_rule (f : W → Y → Z) (hf : Differentiable R (fun (w,y) => f w y)) (g : W → X → Y) (w : W) {g' disc} @@ -247,7 +247,7 @@ theorem comp1_smooth_jumps_rule speed := speed discontinuity := d }) := - comp_smooth_jumps_rule R f g w hf hg + comp_differentiable_jumps_rule R f g w hf hg @[gtrans] @@ -275,7 +275,7 @@ theorem _root_.Prod.mk.arg_fstsnd.HasParamFDerivWithJumpsAt_rule -theorem comp2_smooth_jumps_rule +theorem comp2_differentiable_jumps_rule (f : W → Y₁ → Y₂ → Z) (hf : Differentiable R (fun (w,y₁,y₂) => f w y₁ y₂)) (g₁ : W → X → Y₁) (g₂ : W → X → Y₂) (w : W) {g₁' dg₁} {g₂' dg₂} @@ -303,7 +303,7 @@ theorem comp2_smooth_jumps_rule let y₂ := d.vals x (f w y₁ y₂.1, f w y₁ y₂.2) })) := by - convert comp_smooth_jumps_rule R (fun (w:W) (y:Y₁×Y₂) => f w y.1 y.2) (fun w x => (g₁ w x, g₂ w x)) w + convert comp_differentiable_jumps_rule R (fun (w:W) (y:Y₁×Y₂) => f w y.1 y.2) (fun w x => (g₁ w x, g₂ w x)) w (hf) (Prod.mk.arg_fstsnd.HasParamFDerivWithJumpsAt_rule R g₁ g₂ w hg₁ hg₂ hdisjoint) . simp[Function.comp] @@ -328,49 +328,49 @@ def DisjointJumps {X} [NormedAddCommGroup X] [NormedSpace R X] [MeasureSpace X] @[gtrans] def Prod.fst.arg_self.HasParamFDerivWithJumpsAt_rule := - (comp1_smooth_jumps_rule (R:=R) (W:=W) (X:=X) (Y:=Y×Z) (Z:=Y) (fun _ yz => yz.1) (by fun_prop)) + (comp1_differentiable_jumps_rule (R:=R) (W:=W) (X:=X) (Y:=Y×Z) (Z:=Y) (fun _ yz => yz.1) (by fun_prop)) -- rewrite_type_by (repeat ext); autodiff @[gtrans] def Prod.snd.arg_self.HasParamFDerivWithJumpsAt_rule := - (comp1_smooth_jumps_rule (R:=R) (W:=W) (X:=X) (Y:=Y×Z) (Z:=Z) (fun _ yz => yz.2) (by fun_prop)) + (comp1_differentiable_jumps_rule (R:=R) (W:=W) (X:=X) (Y:=Y×Z) (Z:=Z) (fun _ yz => yz.2) (by fun_prop)) -- rewrite_type_by (repeat ext); autodiff @[gtrans] def HAdd.hAdd.arg_a0a1.HasParamFDerivWithJumpsAt_rule := - (comp2_smooth_jumps_rule (R:=R) (W:=W) (X:=X) (Y₁:=Y) (Y₂:=Y) (Z:=Y) (fun _ y₁ y₂ => y₁ + y₂) (by fun_prop)) + (comp2_differentiable_jumps_rule (R:=R) (W:=W) (X:=X) (Y₁:=Y) (Y₂:=Y) (Z:=Y) (fun _ y₁ y₂ => y₁ + y₂) (by fun_prop)) -- rewrite_type_by (repeat ext); autodiff @[gtrans] def HSub.hSub.arg_a0a1.HasParamFDerivWithJumpsAt_rule := - (comp2_smooth_jumps_rule (R:=R) (W:=W) (X:=X) (Y₁:=Y) (Y₂:=Y) (Z:=Y) (fun _ y₁ y₂ => y₁ - y₂) (by fun_prop)) + (comp2_differentiable_jumps_rule (R:=R) (W:=W) (X:=X) (Y₁:=Y) (Y₂:=Y) (Z:=Y) (fun _ y₁ y₂ => y₁ - y₂) (by fun_prop)) -- rewrite_type_by (repeat ext); autodiff @[gtrans] def Neg.neg.arg_a0.HasParamFDerivWithJumpsAt_rule := - (comp1_smooth_jumps_rule (R:=R) (X:=X) (Y:=Y) (Z:=Y) (fun (w : W) y => - y) (by fun_prop)) + (comp1_differentiable_jumps_rule (R:=R) (X:=X) (Y:=Y) (Z:=Y) (fun (w : W) y => - y) (by fun_prop)) -- rewrite_type_by (repeat ext); autodiff @[gtrans] def HMul.hMul.arg_a0a1.HasParamFDerivWithJumpsAt_rule := - (comp2_smooth_jumps_rule (R:=R) (W:=W) (X:=X) (Y₁:=R) (Y₂:=R) (Z:=R) (fun _ y₁ y₂ => y₁ * y₂) (by fun_prop)) + (comp2_differentiable_jumps_rule (R:=R) (W:=W) (X:=X) (Y₁:=R) (Y₂:=R) (Z:=R) (fun _ y₁ y₂ => y₁ * y₂) (by fun_prop)) -- rewrite_type_by (repeat ext); autodiff @[gtrans] def HPow.hPow.arg_a0.HasParamFDerivWithJumpsAt_rule (n:ℕ) := - (comp1_smooth_jumps_rule (R:=R) (X:=X) (Y:=R) (Z:=R) (fun (w : W) y => y^n) (by fun_prop)) + (comp1_differentiable_jumps_rule (R:=R) (X:=X) (Y:=R) (Z:=R) (fun (w : W) y => y^n) (by fun_prop)) -- rewrite_type_by (repeat ext); autodiff @[gtrans] def HSMul.hSMul.arg_a0a1.HasParamFDerivWithJumpsAt_rule := - (comp2_smooth_jumps_rule (R:=R) (W:=W) (X:=X) (Y₁:=R) (Y₂:=Y) (Z:=Y) (fun _ y₁ y₂ => y₁ • y₂) (by fun_prop)) + (comp2_differentiable_jumps_rule (R:=R) (W:=W) (X:=X) (Y₁:=R) (Y₂:=Y) (Z:=Y) (fun _ y₁ y₂ => y₁ • y₂) (by fun_prop)) -- rewrite_type_by (repeat ext); autodiff @@ -402,7 +402,7 @@ theorem HDiv.hDiv.arg_a0a1.HasParamFDerivWithJumpsAt_rule let z := d.vals x (y/z.1, y/z.2) })) := by - convert comp_smooth_jumps_rule_at (R:=R) + convert comp_differentiable_jumps_rule_at (R:=R) (f:=fun _ (y:R×R) => y.1 / y.2) (g:=fun w x => (f w x, g w x)) (w:=w) (hf:=by simp; sorry_proof) (hg:= Prod.mk.arg_fstsnd.HasParamFDerivWithJumpsAt_rule R f g w hf hg hdisjoint) @@ -422,7 +422,7 @@ theorem ite.arg_te.HasParamFDerivWithJumpsAt_rule (f' := fun dw x => if c w x then f' dw x else g' dw x) (disc := {vals := fun x => (f w x, g w x) - speed := frontierSpeed' R (fun w => {x | c w x}) w + speed := frontierSpeed R (fun w => {x | c w x}) w discontinuity := frontier {x | c w x}} :: df.map (fun d => {d with discontinuity := d.discontinuity ∩ {x | c w x}}) @@ -440,18 +440,18 @@ theorem ite.arg_te.HasParamFDerivWithJumpsAt_rule open Scalar in @[gtrans] def Scalar.sin.arg_a0.HasParamFDerivWithJumpsAt_rule := - (comp1_smooth_jumps_rule (R:=R) (W:=W) (X:=X) (Y:=R) (Z:=R) (fun _ y => sin y) (by simp; fun_prop)) + (comp1_differentiable_jumps_rule (R:=R) (W:=W) (X:=X) (Y:=R) (Z:=R) (fun _ y => sin y) (by simp; fun_prop)) -- rewrite_type_by (repeat ext); autodiff open Scalar in @[gtrans] def Scalar.cos.arg_a0.HasParamFDerivWithJumpsAt_rule := - (comp1_smooth_jumps_rule (R:=R) (W:=W) (X:=X) (Y:=R) (Z:=R) (fun _ y => cos y) (by simp; fun_prop)) + (comp1_differentiable_jumps_rule (R:=R) (W:=W) (X:=X) (Y:=R) (Z:=R) (fun _ y => cos y) (by simp; fun_prop)) -- rewrite_type_by (repeat ext); autodiff @[gtrans] def gaussian.arg_a0.HasParamFDerivWithJumpsAt_rule (σ : R) := - (comp2_smooth_jumps_rule (R:=R) (W:=W) (X:=X) (Y₁:=X) (Y₂:=X) (Z:=R) (fun _ μ x => gaussian μ σ x) (by simp; fun_prop)) + (comp2_differentiable_jumps_rule (R:=R) (W:=W) (X:=X) (Y₁:=X) (Y₂:=X) (Z:=R) (fun _ μ x => gaussian μ σ x) (by simp; fun_prop)) -- rewrite_type_by (repeat ext); autodiff diff --git a/SciLean/Core/Transformations/HasParamDerivWithJumps/HasParamFwdFDerivWithJumps.lean b/SciLean/Core/Transformations/HasParamDerivWithJumps/HasParamFwdFDerivWithJumps.lean index 50edbe94..58e00cde 100644 --- a/SciLean/Core/Transformations/HasParamDerivWithJumps/HasParamFwdFDerivWithJumps.lean +++ b/SciLean/Core/Transformations/HasParamDerivWithJumps/HasParamFwdFDerivWithJumps.lean @@ -108,12 +108,12 @@ theorem smooth_rule unfold HasParamFwdFDerivWithJumpsAt constructor - . apply HasParamFDerivWithJumpsAt.smooth_rule R w f hf + . apply HasParamFDerivWithJumpsAt.differentiable_at_rule R w f hf . simp only [fwdFDeriv, implies_true] -theorem comp_smooth_jumps_rule +theorem comp_differentiable_jumps_rule (f : W → Y → Z) (g : W → X → Y) (w : W) {g' disc} (hf : Differentiable R (fun (w,y) => f w y)) @@ -133,7 +133,7 @@ theorem comp_smooth_jumps_rule unfold HasParamFwdFDerivWithJumpsAt constructor - . convert HasParamFDerivWithJumpsAt.comp_smooth_jumps_rule R f g w hf hg.1 + . convert HasParamFDerivWithJumpsAt.comp_differentiable_jumps_rule R f g w hf hg.1 simp [fwdFDeriv, hg.2] . simp [fwdFDeriv,hg.2] @@ -145,8 +145,7 @@ theorem _root_.Prod.mk.arg_fstsnd.HasParamFwdFDerivWithJumpsAt_rule {f' fdisc} {g' gdisc} (hf : HasParamFwdFDerivWithJumpsAt R f w f' fdisc) (hg : HasParamFwdFDerivWithJumpsAt R g w g' gdisc) - (hdisjoint : AlmostDisjoint fdisc.getDiscontinuity gdisc.getDiscontinuity μH[finrank ℝ X - (1:ℕ)]) - /- (hIJ : DisjointJumps R Sf Sg) -/ : + (hdisjoint : AlmostDisjoint fdisc.getDiscontinuity gdisc.getDiscontinuity μH[finrank ℝ X - (1:ℕ)]) : HasParamFwdFDerivWithJumpsAt (R:=R) (fun w x => (f w x, g w x)) w (f' := fun dw x => let ydy := f' dw x @@ -171,7 +170,7 @@ theorem _root_.Prod.mk.arg_fstsnd.HasParamFwdFDerivWithJumpsAt_rule . simp [hf.2, hg.2] -theorem comp1_smooth_jumps_rule +theorem comp1_differentiable_jumps_rule (f : W → Y → Z) (hf : Differentiable R (fun (w,y) => f w y)) (g : W → X → Y) (w : W) {g' disc} @@ -188,10 +187,10 @@ theorem comp1_smooth_jumps_rule speed := speed discontinuity := d }) := - comp_smooth_jumps_rule f g w hf hg + comp_differentiable_jumps_rule f g w hf hg -theorem comp2_smooth_jumps_rule +theorem comp2_differentiable_jumps_rule (f : W → Y₁ → Y₂ → Z) (hf : Differentiable R (fun (w,y₁,y₂) => f w y₁ y₂)) (g₁ : W → X → Y₁) (g₂ : W → X → Y₂) (w : W) {g₁' dg₁} {g₂' dg₂} @@ -217,7 +216,7 @@ theorem comp2_smooth_jumps_rule let y₂ := d.vals x (f w y₁ y₂.1, f w y₁ y₂.2) })) := by - convert comp_smooth_jumps_rule (R:=R) (fun w (y:Y₁×Y₂) => f w y.1 y.2) (fun w x => (g₁ w x, g₂ w x)) w + convert comp_differentiable_jumps_rule (R:=R) (fun w (y:Y₁×Y₂) => f w y.1 y.2) (fun w x => (g₁ w x, g₂ w x)) w hf (by gtrans (disch:=first | fun_prop | assumption)) . simp[List.map_append]; rfl @@ -228,42 +227,42 @@ open HasParamFwdDerivWithJumps @[gtrans] def Prod.fst.arg_self.HasParamFwdFDerivWithJumpsAt_rule := - (comp1_smooth_jumps_rule (R:=R) (W:=W) (X:=X) (Y:=Y×Z) (Z:=Y) (fun _ yz => yz.1) (by fun_prop)) + (comp1_differentiable_jumps_rule (R:=R) (W:=W) (X:=X) (Y:=Y×Z) (Z:=Y) (fun _ yz => yz.1) (by fun_prop)) rewrite_type_by (repeat ext); autodiff @[gtrans] def Prod.snd.arg_self.HasParamFwdFDerivWithJumpsAt_rule := - (comp1_smooth_jumps_rule (R:=R) (W:=W) (X:=X) (Y:=Y×Z) (Z:=Z) (fun _ yz => yz.2) (by fun_prop)) + (comp1_differentiable_jumps_rule (R:=R) (W:=W) (X:=X) (Y:=Y×Z) (Z:=Z) (fun _ yz => yz.2) (by fun_prop)) rewrite_type_by (repeat ext); autodiff @[gtrans] def HAdd.hAdd.arg_a0a1.HasParamFwdFDerivWithJumpsAt_rule := - (comp2_smooth_jumps_rule (R:=R) (W:=W) (X:=X) (Y₁:=Y) (Y₂:=Y) (Z:=Y) (fun _ y₁ y₂ => y₁ + y₂) (by fun_prop)) + (comp2_differentiable_jumps_rule (R:=R) (W:=W) (X:=X) (Y₁:=Y) (Y₂:=Y) (Z:=Y) (fun _ y₁ y₂ => y₁ + y₂) (by fun_prop)) rewrite_type_by (repeat ext); autodiff @[gtrans] def HSub.hSub.arg_a0a1.HasParamFwdFDerivWithJumpsAt_rule := - (comp2_smooth_jumps_rule (R:=R) (W:=W) (X:=X) (Y₁:=Y) (Y₂:=Y) (Z:=Y) (fun _ y₁ y₂ => y₁ - y₂) (by fun_prop)) + (comp2_differentiable_jumps_rule (R:=R) (W:=W) (X:=X) (Y₁:=Y) (Y₂:=Y) (Z:=Y) (fun _ y₁ y₂ => y₁ - y₂) (by fun_prop)) rewrite_type_by (repeat ext); autodiff @[gtrans] def Neg.neg.arg_a0.HasParamFwdFDerivWithJumpsAt_rule := - (comp1_smooth_jumps_rule (R:=R) (X:=X) (Y:=Y) (Z:=Y) (fun (w : W) y => - y) (by fun_prop)) + (comp1_differentiable_jumps_rule (R:=R) (X:=X) (Y:=Y) (Z:=Y) (fun (w : W) y => - y) (by fun_prop)) rewrite_type_by (repeat ext); autodiff @[gtrans] def HMul.hMul.arg_a0a1.HasParamFwdFDerivWithJumpsAt_rule := - (comp2_smooth_jumps_rule (R:=R) (W:=W) (X:=X) (Y₁:=R) (Y₂:=R) (Z:=R) (fun _ y₁ y₂ => y₁ * y₂) (by fun_prop)) + (comp2_differentiable_jumps_rule (R:=R) (W:=W) (X:=X) (Y₁:=R) (Y₂:=R) (Z:=R) (fun _ y₁ y₂ => y₁ * y₂) (by fun_prop)) rewrite_type_by (repeat ext); autodiff @[gtrans] def HPow.hPow.arg_a0.HasParamFwdFDerivWithJumpsAt_rule (n:ℕ) := - (comp1_smooth_jumps_rule (R:=R) (X:=X) (Y:=R) (Z:=R) (fun (w : W) y => y^n) (by fun_prop)) + (comp1_differentiable_jumps_rule (R:=R) (X:=X) (Y:=R) (Z:=R) (fun (w : W) y => y^n) (by fun_prop)) rewrite_type_by (repeat ext); autodiff @[gtrans] def HSMul.hSMul.arg_a0a1.HasParamFwdFDerivWithJumpsAt_rule := - (comp2_smooth_jumps_rule (R:=R) (W:=W) (X:=X) (Y₁:=R) (Y₂:=Y) (Z:=Y) (fun _ y₁ y₂ => y₁ • y₂) (by fun_prop)) + (comp2_differentiable_jumps_rule (R:=R) (W:=W) (X:=X) (Y₁:=R) (Y₂:=Y) (Z:=Y) (fun _ y₁ y₂ => y₁ • y₂) (by fun_prop)) rewrite_type_by (repeat ext); autodiff @[gtrans] @@ -313,7 +312,7 @@ theorem ite.arg_te.HasParamFwdFDerivWithJumpsAt_rule (f' := fun dw x => if c w x then f' dw x else g' dw x) (disc := {vals := fun x => (f w x, g w x) - speed := frontierSpeed' R (fun w => {x | c w x}) w + speed := frontierSpeed R (fun w => {x | c w x}) w discontinuity := frontier {x | c w x}} :: df.map (fun d => {d with discontinuity := d.discontinuity ∩ {x | c w x}}) @@ -337,18 +336,18 @@ theorem ite.arg_te.HasParamFwdFDerivWithJumpsAt_rule open Scalar in @[gtrans] def Scalar.sin.arg_a0.HasParamFwdFDerivWithJumpsAt_rule := - (comp1_smooth_jumps_rule (R:=R) (W:=W) (X:=X) (Y:=R) (Z:=R) (fun _ y => sin y) (by simp; fun_prop)) + (comp1_differentiable_jumps_rule (R:=R) (W:=W) (X:=X) (Y:=R) (Z:=R) (fun _ y => sin y) (by simp; fun_prop)) rewrite_type_by (repeat ext); autodiff open Scalar in @[gtrans] def Scalar.cos.arg_a0.HasParamFwdFDerivWithJumpsAt_rule := - (comp1_smooth_jumps_rule (R:=R) (W:=W) (X:=X) (Y:=R) (Z:=R) (fun _ y => cos y) (by simp; fun_prop)) + (comp1_differentiable_jumps_rule (R:=R) (W:=W) (X:=X) (Y:=R) (Z:=R) (fun _ y => cos y) (by simp; fun_prop)) rewrite_type_by (repeat ext); autodiff @[gtrans] def gaussian.arg_a0.HasParamFwdFDerivWithJumpsAt_rule (σ : R) := - (comp2_smooth_jumps_rule (R:=R) (W:=W) (X:=X) (Y₁:=X) (Y₂:=X) (Z:=R) (fun _ μ x => gaussian μ σ x) (by simp; fun_prop)) + (comp2_differentiable_jumps_rule (R:=R) (W:=W) (X:=X) (Y₁:=X) (Y₂:=X) (Z:=R) (fun _ μ x => gaussian μ σ x) (by simp; fun_prop)) rewrite_type_by (repeat ext); autodiff diff --git a/SciLean/Core/Transformations/HasParamDerivWithJumps/HasParamRevFDerivWithJumps.lean b/SciLean/Core/Transformations/HasParamDerivWithJumps/HasParamRevFDerivWithJumps.lean index 53f87eb7..5b4bbf41 100644 --- a/SciLean/Core/Transformations/HasParamDerivWithJumps/HasParamRevFDerivWithJumps.lean +++ b/SciLean/Core/Transformations/HasParamDerivWithJumps/HasParamRevFDerivWithJumps.lean @@ -24,11 +24,13 @@ set_default_scalar R ---------------------------------------------------------------------------------------------------- +@[fun_trans] theorem adjoint_integral (f : X → Y → Z) (μ : Measure X) : adjoint R (fun y => ∫ x, f x y ∂μ) = fun y' => ∫ x, adjoint R (f x) y' ∂μ := sorry_proof +@[fun_trans] theorem adjoint_sum {ι} [IndexType ι] (f : ι → Y → Z) : adjoint R (fun y => ∑ i, f i y) = @@ -62,8 +64,8 @@ attribute [fun_trans] adjoint_adjoint variable (R) open Classical in noncomputable -def frontierGrad (A : W → Set X) (w : W) (x : X) : W := - adjoint R (fun dw => frontierSpeed' R A w dw x) 1 +def frontierGrad (Ω : W → Set X) (w : W) (x : X) : W := + adjoint R (fun dw => frontierSpeed R Ω w dw x) 1 variable (W X Y) @@ -88,7 +90,7 @@ def HasParamRevFDerivWithJumpsAt (f : W → X → Y) (w : W) (f' : outParam <| X → Y×(Y→W)) (disc : outParam <| DiscontinuityRevDataList W X Y) := HasParamFDerivWithJumpsAt R f w - (fun (dw : W) (x : X) => adjoint R (fun dy => ⟪(f' x).2 dy, dw⟫) 1) + (fun (dw : W) (x : X) => adjoint R (f' x).2 dw) (disc.map (fun ⟨df,s,S⟩ => ⟨df,fun dy x => ⟪s x, dy⟫,S⟩)) ∧ ∀ x, (f' x).1 = f w x @@ -126,25 +128,25 @@ theorem revFDeriv_under_integral @[fun_trans] theorem revFDeriv_under_integral_over_set - (f : W → X → Y) (w : W) (μ : Measure X) (A : Set X) + (f : W → X → Y) (w : W) (μ : Measure X) (Ω : Set X) {f' disc} (hf : HasParamRevFDerivWithJumpsAt R f w f' disc) - (hA : AlmostDisjoint (frontier A) disc.getDiscontinuity μH[finrank ℝ X - (1:ℕ)]) : - (revFDeriv R (fun w' => ∫ x in A, f w' x ∂μ) w) + (hA : AlmostDisjoint (frontier Ω) disc.getDiscontinuity μH[finrank ℝ X - (1:ℕ)]) : + (revFDeriv R (fun w' => ∫ x in Ω, f w' x ∂μ) w) = - let val := ∫ x in A, f w x ∂μ + let val := ∫ x in Ω, f w x ∂μ (val, fun dy => - let interior := ∫ x in A, (f' x).2 dy ∂μ + let interior := ∫ x in Ω, (f' x).2 dy ∂μ let density := fun x => Scalar.ofENNReal (R:=R) (μ.rnDeriv volume x) let shocks := disc.foldl (init:=0) - fun sum ⟨df,s,S⟩ => sum + - ∫ x in S ∩ A, + fun sum ⟨df,s,Γ⟩ => sum + + ∫ x in Γ ∩ Ω, let vals := df x (⟪vals.1 - vals.2, dy⟫ * density x) • s x ∂μH[finrank R X - (1:ℕ)] interior + shocks) := by simp[revFDeriv] - simp only [fderiv_under_integral_over_set R f w _ μ A hf.1 sorry_proof] + simp only [fderiv_under_integral_over_set R f w _ μ Ω hf.1 sorry_proof] have hf' : ∀ x, IsContinuousLinearMap R (f' x).2 := sorry_proof -- this should be part of hf fun_trans (disch:=apply hf') [adjoint_sum,adjoint_integral,adjoint_adjoint,smul_smul] @@ -166,12 +168,12 @@ theorem smooth_rule unfold HasParamRevFDerivWithJumpsAt constructor - . convert HasParamFDerivWithJumpsAt.smooth_rule R w f hf + . convert HasParamFDerivWithJumpsAt.differentiable_at_rule R w f hf . fun_trans [revFDeriv] . simp [revFDeriv] -theorem comp_smooth_jumps_rule +theorem comp_differentiable_jumps_rule (f : W → Y → Z) (g : W → X → Y) (w : W) {g' disc} (hf : Differentiable R (fun (w,y) => f w y)) @@ -194,7 +196,7 @@ theorem comp_smooth_jumps_rule unfold HasParamRevFDerivWithJumpsAt constructor - . convert HasParamFDerivWithJumpsAt.comp_smooth_jumps_rule R f g w hf hg.1 + . convert HasParamFDerivWithJumpsAt.comp_differentiable_jumps_rule R f g w hf hg.1 . rename_i w x have hg' : IsContinuousLinearMap R (g' x).2 := by sorry_proof simp [revFDeriv, hg.2] @@ -243,7 +245,7 @@ theorem _root_.Prod.mk.arg_fstsnd.HasParamRevFDerivWithJumpsAt_rule -theorem comp1_smooth_jumps_rule +theorem comp1_differentiable_jumps_rule (f : W → Y → Z) (hf : Differentiable R (fun (w,y) => f w y)) (g : W → X → Y) (w : W) {g' disc} @@ -264,11 +266,11 @@ theorem comp1_smooth_jumps_rule speedGrad := speedGrad discontinuity := d }) := - comp_smooth_jumps_rule f g w hf hg + comp_differentiable_jumps_rule f g w hf hg -theorem comp2_smooth_jumps_rule +theorem comp2_differentiable_jumps_rule (f : W → Y₁ → Y₂ → Z) (hf : Differentiable R (fun (w,y₁,y₂) => f w y₁ y₂)) (g₁ : W → X → Y₁) (g₂ : W → X → Y₂) (w : W) {g₁' dg₁} {g₂' dg₂} @@ -298,7 +300,7 @@ theorem comp2_smooth_jumps_rule let y₂ := d.vals x (f w y₁ y₂.1, f w y₁ y₂.2) })) := by - convert comp_smooth_jumps_rule (R:=R) (fun w (y:Y₁×Y₂) => f w y.1 y.2) (fun w x => (g₁ w x, g₂ w x)) w + convert comp_differentiable_jumps_rule (R:=R) (fun w (y:Y₁×Y₂) => f w y.1 y.2) (fun w x => (g₁ w x, g₂ w x)) w hf (by gtrans (disch:=first | fun_prop | assumption)) . fun_trans [hg₁.2,hg₂.2]; ac_rfl . simp[List.map_append]; rfl @@ -310,42 +312,42 @@ open HasParamRevFDerivWithJumpsAt @[gtrans] def Prod.fst.arg_self.HasParamRevFDerivWithJumpsAt_rule := - (comp1_smooth_jumps_rule (R:=R) (W:=W) (X:=X) (Y:=Y×Z) (Z:=Y) (fun _ yz => yz.1) (by fun_prop)) + (comp1_differentiable_jumps_rule (R:=R) (W:=W) (X:=X) (Y:=Y×Z) (Z:=Y) (fun _ yz => yz.1) (by fun_prop)) rewrite_type_by (repeat ext); autodiff @[gtrans] def Prod.snd.arg_self.HasParamRevFDerivWithJumpsAt_rule := - (comp1_smooth_jumps_rule (R:=R) (W:=W) (X:=X) (Y:=Y×Z) (Z:=Z) (fun _ yz => yz.2) (by fun_prop)) + (comp1_differentiable_jumps_rule (R:=R) (W:=W) (X:=X) (Y:=Y×Z) (Z:=Z) (fun _ yz => yz.2) (by fun_prop)) rewrite_type_by (repeat ext); autodiff @[gtrans] def HAdd.hAdd.arg_a0a1.HasParamRevFDerivWithJumpsAt_rule := - (comp2_smooth_jumps_rule (R:=R) (W:=W) (X:=X) (Y₁:=Y) (Y₂:=Y) (Z:=Y) (fun _ y₁ y₂ => y₁ + y₂) (by fun_prop)) + (comp2_differentiable_jumps_rule (R:=R) (W:=W) (X:=X) (Y₁:=Y) (Y₂:=Y) (Z:=Y) (fun _ y₁ y₂ => y₁ + y₂) (by fun_prop)) rewrite_type_by (repeat ext); autodiff @[gtrans] def HSub.hSub.arg_a0a1.HasParamRevFDerivWithJumpsAt_rule := - (comp2_smooth_jumps_rule (R:=R) (W:=W) (X:=X) (Y₁:=Y) (Y₂:=Y) (Z:=Y) (fun _ y₁ y₂ => y₁ - y₂) (by fun_prop)) + (comp2_differentiable_jumps_rule (R:=R) (W:=W) (X:=X) (Y₁:=Y) (Y₂:=Y) (Z:=Y) (fun _ y₁ y₂ => y₁ - y₂) (by fun_prop)) rewrite_type_by (repeat ext); autodiff @[gtrans] def Neg.neg.arg_a0.HasParamRevFDerivWithJumpsAt_rule := - (comp1_smooth_jumps_rule (R:=R) (W:=W) (X:=X) (Y:=Y) (Z:=Y) (fun _ y => - y) (by fun_prop)) + (comp1_differentiable_jumps_rule (R:=R) (W:=W) (X:=X) (Y:=Y) (Z:=Y) (fun _ y => - y) (by fun_prop)) rewrite_type_by (repeat ext); autodiff @[gtrans] def HMul.hMul.arg_a0a1.HasParamRevFDerivWithJumpsAt_rule := - (comp2_smooth_jumps_rule (R:=R) (W:=W) (X:=X) (Y₁:=R) (Y₂:=R) (Z:=R) (fun _ y₁ y₂ => y₁ * y₂) (by fun_prop)) + (comp2_differentiable_jumps_rule (R:=R) (W:=W) (X:=X) (Y₁:=R) (Y₂:=R) (Z:=R) (fun _ y₁ y₂ => y₁ * y₂) (by fun_prop)) rewrite_type_by (repeat ext); autodiff @[gtrans] def HPow.hPow.arg_a0.HasParamRevFDerivWithJumpsAt_rule (n:ℕ) := - (comp1_smooth_jumps_rule (R:=R) (X:=X) (Y:=R) (Z:=R) (fun (w : W) y => y^n) (by fun_prop)) + (comp1_differentiable_jumps_rule (R:=R) (X:=X) (Y:=R) (Z:=R) (fun (w : W) y => y^n) (by fun_prop)) rewrite_type_by (repeat ext); autodiff @[gtrans] def HSMul.hSMul.arg_a0a1.HasParamRevFDerivWithJumpsAt_rule := - (comp2_smooth_jumps_rule (R:=R) (W:=W) (X:=X) (Y₁:=R) (Y₂:=Y) (Z:=Y) (fun _ y₁ y₂ => y₁ • y₂) (by fun_prop)) + (comp2_differentiable_jumps_rule (R:=R) (W:=W) (X:=X) (Y₁:=R) (Y₂:=Y) (Z:=Y) (fun _ y₁ y₂ => y₁ • y₂) (by fun_prop)) rewrite_type_by (repeat ext); autodiff @[gtrans] @@ -426,18 +428,18 @@ theorem ite.arg_te.HasParamRevFDerivWithJumpsAt_rule open Scalar in @[gtrans] def Scalar.sin.arg_a0.HasParamRevFDerivWithJumpsAt_rule := - (comp1_smooth_jumps_rule (R:=R) (W:=W) (X:=X) (Y:=R) (Z:=R) (fun _ y => sin y) (by simp; fun_prop)) + (comp1_differentiable_jumps_rule (R:=R) (W:=W) (X:=X) (Y:=R) (Z:=R) (fun _ y => sin y) (by simp; fun_prop)) rewrite_type_by (repeat ext); autodiff open Scalar in @[gtrans] def Scalar.cos.arg_a0.HasParamRevFDerivWithJumpsAt_rule := - (comp1_smooth_jumps_rule (R:=R) (W:=W) (X:=X) (Y:=R) (Z:=R) (fun _ y => cos y) (by simp; fun_prop)) + (comp1_differentiable_jumps_rule (R:=R) (W:=W) (X:=X) (Y:=R) (Z:=R) (fun _ y => cos y) (by simp; fun_prop)) rewrite_type_by (repeat ext); autodiff @[gtrans] def gaussian.arg_a0.HasParamRevFDerivWithJumpsAt_rule (σ : R) := - (comp2_smooth_jumps_rule (R:=R) (W:=W) (X:=X) (Y₁:=X) (Y₂:=X) (Z:=R) (fun _ μ x => gaussian μ σ x) (by simp; fun_prop)) + (comp2_differentiable_jumps_rule (R:=R) (W:=W) (X:=X) (Y₁:=X) (Y₂:=X) (Z:=R) (fun _ μ x => gaussian μ σ x) (by simp; fun_prop)) rewrite_type_by (repeat ext); autodiff diff --git a/SciLean/Doodle.lean b/SciLean/Doodle.lean index c3b13c8d..9532342e 100644 --- a/SciLean/Doodle.lean +++ b/SciLean/Doodle.lean @@ -1,216 +1,42 @@ -import SciLean.Tactic.IfPull - +import SciLean open SciLean -set_option trace.Meta.Tactic.if_pull true - -#check (fun x : Float => if 1 ≤ 0 then if x ≤ 1 then 1 else 0 else 0) rewrite_by - simp only [Tactic.if_pull] - - - -#check (let y := 5; (if 0 < 1 then (fun x : Float => x + 2 + y) else (fun x : Float => x + 3 + y)) 42).log rewrite_by - simp (config:={zeta:=false}) only [Tactic.if_pull] - -#exit - -open SciLean - -set_default_scalar Float - - - -#check ∂>! x : Float, ((x*x)*x) - -#check HMul.hMul.arg_a0a1. - - -#check ∂! x : Float, x*x - -def foo (x : Float) := ∂! (x':=x), x'*x' - -#exit - -def dot {n : Nat} (x y : Float^[n]) : Float := ∑ i, x[i] * y[i] - -variable (a b : Int) -#check DataArrayN Float (Set.Icc a b) - -open Set - -#check Icc (-1) 1 - -def ab := ⊞ (i : Fin 10) => 1.0 -def _root_.Int.toFloat (a : Int) : Float := - if a ≥ 0 then - a.toNat.toFloat - else - -(-a).toNat.toFloat - -#check (↑(Icc (-1) 1) : Type) - -#check Elem --- set_option pp.universes true in --- set_option pp.coercions false in --- set_option pp.all true in - -def D := fun (i j : Icc (-1) 1) => (i.1.toFloat + j.1.toFloat) - -#synth ArrayTypeNotation _ ((Icc (-1) 5)×Fin 10) Float -#synth IndexType ((Icc (-1) 2) × Fin 10) -#synth ArrayTypeNotation _ ((Icc (-1) 2) × Fin 10) Float - -#check ⊞ (i j : Fin 10) => (1.0 : Float) -#check ⊞ (i j : Icc (-1) 1) => (1.0 : Float) - -#eval dot ⊞[1.0,1.0] ⊞[1.0,1.0] -#eval dot ⊞[1.0,1.0] ⊞[1.0,1.0,1.0] - -def u := ⊞[(1.0 : Float), 2.0] - -#check u - -#eval u[0] -#eval u[1] -#eval fun i j => (u[i],u[j]) - -#check ⊞[⊞[(1.0:Float), 2.0],⊞[(3.0:Float), 4.0]] - -× -def A := ⊞[1.0, 2.0; - 3.0, 4.0] - - -instance : Coe Nat Float := ⟨fun i => i.toUInt64.toFloat⟩ -def _root_.Fin.toFloat {n} (i : Fin n) : Float := i.1.toUInt64.toFloat - -#check ⊞ (i : Fin 10) => i.toFloat - -variable {n m: Nat} (x : Float^[n]) - -#check ⊞ i (j : Fin m) => x[i]^j.1 - -variable (A : Float^[(2 : Nat),(2 : Nat)]) - --- #check A[0,1] --- #check Float^[[2,2], [2,2]] -- Fin^[4,4] - -def outerProduct {n m : Nat} (x : Float^[n]) (y : Float^[m]) : Float^[n,m] := - ⊞ i j => x[i]*y[j] - -def outerProduct' {n m : Nat} (x : Float^[n]) (y : Float^[m]) : (Float^[m])^[n] := - ⊞ i => (⊞ j => x[i]*y[j]) - ∑ j, A[i,j] * x[j] - - -variable (n : Nat) -#synth VAdd ℤ (Fin n) - -instance (n : Nat) : VAdd ℤ (Fin n) := ⟨i.1 +ᵥ j⟩ - - -def shift ( - - -#check Float^[[-1:1]] - - -#check Fold.fold - -#check ArrayType.max - -#eval 1.0 ⊔ 2.0 -#check Lattice - - -open Scalar -def softMax {I} [IndexType I] - (r : Float) (x : Float^I) : Float^I := Id.run do - let m := x.max - - let mut w := 0.0 - let mut x := x - - for i in IndexType.univ I do - let xi := exp (r*(x[i]-m)) - x[i] := xi - w += xi - - x := (1/w) • x - return x - -- let x := ArrayType.map (fun xi => Scalar.exp (r*(xi-m))) x - -- let w := ∑ i, x[i] - -- (1/w) • x - -#check ArrayType.mapMono - - - -open Scalar - -#eval DataArrayN.mapMono (⊞ (i j k : Fin 2) => (IndexType.toFin (i,j,k)).toFloat) (fun x => sqrt x) - - - -#eval (0 : Float^[3]) |>.mapIdxMono (fun i _ => i.toFloat) |>.mapMono (fun x => sqrt x) - -#eval ⊞[(1.0 :Float),2.0,3.0].fold (·+·) 0 - - --- n#eval (·+·) -- - -def x' := ⊞[(1.0 : Float),2.0,3.0,4.0] -def A' := ⊞ (i j : Fin 2) => (1 + 2*i.toFloat + 4*j.toFloat) - -#eval ∑ i, x'[i] -#eval ∏ i, x'[i] - -#eval ∑ i j, A'[i,j] eval ∏ i j, A'[i,j] - - - - - - - -def _root_.Fin.shift {n} (i : Fin n) (j : ℤ) : Fin n := - { val := ((i.1 + j) % n ).toNat, isLt := sorry } - - - - -def conv1d {n k : Nat} (x : Float^[n]) (w : Float^[[-k:k]]) := - ⊞ (i : Fin n) => ∑ j, w[j] * x[i.shift j.1] - - -def conv2d {n m k : Nat} (x : Float^[n,m]) (w : Float^[[-k:k],[-k:k]]) := - ⊞ (i : Fin n) (j : Fin m) => ∑ a b, w[a,b] * x[i.shift a, j.shift b] +def shift {X} [Add X] (u : X) (f : X → Y) (x : X) := f (x + u) +variable + {W} [NormedAddCommGroup W] [NormedSpace ℝ W] -def conv2d' {n m : Nat} (k : Nat) (J : Type) {I : Type} [IndexType I] [IndexType J] [DecidableEq J] - (w : Float^[I,J,[-k:k],[-k:k]]) (b : Float^[J,n,m]) (x : Float^[I,n,m]) : Float^[J,n,m] := - ⊞ κ (i : Fin n) (j : Fin m) => ∑ ι a b, w[ι,κ,a,b] * x[ι, i.shift a, j.shift b] +set_option trace.Meta.Tactic.fun_trans true in +set_option trace.Meta.Tactic.simp.rewrite true in +theorem shift_fwdFDeriv (t : W → ℝ) (f : W → ℝ → ℝ) (x : W → ℝ) + (ht : Differentiable ℝ t) (hf : Differentiable ℝ (fun (w,r) => f w r)) + (hx : Differentiable ℝ x) : + fwdFDeriv ℝ (fun w => shift (t w) (f w) (x w)) + = + fun w dw => + let df := fun x dx => fwdFDeriv ℝ (fun (w',x') => f w' x') (w,x) (dw,dx) + let (tw,dt) := fwdFDeriv ℝ t w dw + let (xw,dx) := fwdFDeriv ℝ x w dw + shift tw (df · (dx + dt)) xw := by + funext w dw + unfold shift + conv => lhs; lautodiff + simp -variable (x : Float^[3, 10,10]) +#check Prod.fst_add -#check - fun (w₁,b₁,w₂,b₂) => - x |> conv2d' 3 (Fin 3) w₁ b₁ - |> conv2d' 5 (Fin 10) w₂ b₂ -#check fun w b => conv2d' 3 (Fin 3) (I:=Fin 2) w +#print axioms shift_fwdFDeriv -#check conv2d' -instance : VAdd ℤ (Fin n) := ⟨fun j i => i.shift j⟩ +#print axioms SciLean.FwdFDeriv.apply_rule -def conv1d' {n k : Nat} (x : Float^[n]) (w : Float^[[-k:k]]) := - ⊞ (i : Fin n) => ∑ j, w[j] * x[j.1 +ᵥ i] +#print axioms SciLean.FwdFDeriv.Prod.mk.arg_fstsnd.fwdFDeriv_rule_at -#check Float^[3, 4] -#check Float^[Fin 3, Fin 4] +#print axioms SciLean.FwdFDeriv.comp_rule_at -#eval (-14 : ℤ) % 10 +#print axioms SciLean.FwdFDeriv.HAdd.hAdd.arg_a0a1.fwdFDeriv_rule_at +#print axioms SciLean.FwdFDeriv.id_rule diff --git a/test/param_deriv/affine_discont_2d.lean b/test/param_deriv/affine_discont_2d.lean index dc944253..1cd2db74 100644 --- a/test/param_deriv/affine_discont_2d.lean +++ b/test/param_deriv/affine_discont_2d.lean @@ -75,5 +75,5 @@ def test_fgradient (numSamples : ℕ) (a b c d : R) (w : R) := #eval 0 -#eval (test_fderiv 1000 1.0 1.0 0.0 1.0 0.3).get -#eval (test_fgradient 1000 1.0 1.0 0.0 1.0 0.3).get +#eval (test_fderiv 1000 1.0 1.0 0.3 1.0 0.5).get +-- #eval (test_fgradient 1000 1.0 1.0 0.0 1.0 0.3).get diff --git a/test/param_deriv/circle_discont_2d.lean b/test/param_deriv/circle_discont_2d.lean index 43c2faf2..376aed18 100644 --- a/test/param_deriv/circle_discont_2d.lean +++ b/test/param_deriv/circle_discont_2d.lean @@ -27,7 +27,7 @@ theorem circle_parametric_inverse_polar (r : R) (hr : r ≠ 0) : SurfaceParametrization {x : R×R | x.1 ^ 2 + x.2 ^ 2 = r ^ 2} (U := R) (param := - [(Ico 0 (2*π), circleParam)]) := sorry_proof + [(Ico 0 (2*π), fun x => r•circleParam x)]) := sorry_proof def test_fderiv (numSamples : ℕ) (w : R) := @@ -68,6 +68,4 @@ def test_fderiv (numSamples : ℕ) (w : R) := - -#eval (test_fderiv 1000 (-1.5)).get -#eval (test_fderiv 1000 (1.5)).get +#eval Rand.print_mean_variance (test_fderiv 1 0.5) 5000 "" diff --git a/test/param_deriv/risk_minimization.lean b/test/param_deriv/risk_minimization.lean index 969b2025..762cb19f 100644 --- a/test/param_deriv/risk_minimization.lean +++ b/test/param_deriv/risk_minimization.lean @@ -27,11 +27,7 @@ def test_fwdFDeriv (l u : R) (θ dθ : R×R×R) := (norm:=lsimp only [ftrans_simp])) (hA := by assume_almost_disjoint)] - lautodiff (disch:=first | fun_prop | gtrans (disch:=fun_prop)) - - conv in (∫ _ in _, _ ∂μH[0]) => - lautodiff (disch:=first | assumption | fun_prop | apply hθ) - + lautodiff (disch:=first | fun_prop | gtrans (disch:=fun_prop) | assumption) conv in (occs:=*) (∫ _ in _, _ ∂_) => . lsimp only [Rand.integral_eq_uniform_E R, @@ -41,11 +37,30 @@ def test_fwdFDeriv (l u : R) (θ dθ : R×R×R) := pull_mean -#exit +-- def test_fgradient (l u : R) (θ : R×R×R) := +-- derive_random_approx +-- (fgradient (fun ((a,b,μ) : R×R×R) => ∫ x in Icc l u, +-- ‖ if x ∈ Icc a b then gaussian μ (5:R) x else 0 - gaussian 2 (5:R) x ‖₂²) θ) +-- assuming (hθ : θ.1 < θ.2.1) +-- by +-- unfold fgradient +-- rw[revFDeriv_under_integral_over_set +-- (hf:= by gtrans +-- (disch:=first | fun_prop_no_ifs | assume_almost_disjoint) +-- (norm:=lsimp only [ftrans_simp])) +-- (hA := by assume_almost_disjoint)] + +-- lautodiff (disch:=first | fun_prop | gtrans (disch:=fun_prop) | assumption) + +-- conv in (occs:=*) (∫ _ in _, _ ∂_) => +-- . lsimp only [Rand.integral_eq_uniform_E R, +-- Rand.E_eq_mean_estimateE R 10] +-- lsimp only [ftrans_simp] + +-- pull_mean + #eval 0 -#eval test_fderiv 0.1 -#eval (test_fwdFDeriv 100 0.1).get -#eval test_fgradient 0.1 +#eval (test_fwdFDeriv (-1.0) 1.0 (0.0,1.0,0.5) (1.0,0.0,0.0)).get diff --git a/test/param_deriv/simple_discont_1d.lean b/test/param_deriv/simple_discont_1d.lean index 09d583e8..0628e4ac 100644 --- a/test/param_deriv/simple_discont_1d.lean +++ b/test/param_deriv/simple_discont_1d.lean @@ -64,6 +64,6 @@ def test_fgradient (w : R) := #eval 0 -#eval test_fderiv 0.1 +#eval test_fderiv 0.5 #eval (test_fwdFDeriv 100 0.1).get #eval test_fgradient 0.1 diff --git a/test/param_deriv/simple_discont_2d.lean b/test/param_deriv/simple_discont_2d.lean index f73857a7..7b1055ac 100644 --- a/test/param_deriv/simple_discont_2d.lean +++ b/test/param_deriv/simple_discont_2d.lean @@ -18,7 +18,7 @@ def test_fderiv (numSamples : ℕ) (w : R) := derive_random_approx (fderiv R (fun w' => ∫ xy in Icc (0:R) 1 ×ˢ (Icc (0 : R) 1), - if xy.1 ≤ w' then (1:R) else (0:R)) w 1) + if xy.1 + xy.2 ≤ w' then (1:R) else (0:R)) w 1) by rw[fderiv_under_integral_over_set (hf:= by gtrans @@ -103,6 +103,6 @@ def test_fgrad (numSamples : ℕ) (w : R) := #eval 0 -#eval (test_fderiv 100 0.5).get +#eval Rand.print_mean_variance (test_fderiv 1 0.5) 1000 "" #eval (test_fwdFDeriv 100 0.5).get #eval (test_fgrad 100 0.5).get