/SphericalCoord_aristotle.lean at c314c45fd1c615c3e67faacc0edab1a21d5fb5bb

/-
This file was edited by Aristotle.

Lean Toolchain version: leanprover/lean4:v4.20.0-rc5
Mathlib version: d62eab0cc36ea522904895389c301cf8d844fd69 (May 9, 2025)

Your Lean code is run in a custom environment, which uses these headers:

set_option maxHeartbeats 0
set_option maxRecDepth 4000
set_option synthInstance.maxHeartbeats 20000
set_option synthInstance.maxSize 128

The following was proved by Aristotle:

- def sphericalCoord : PartialHomeomorph (ℝ × ℝ × ℝ) (ℝ × ℝ × ℝ) where
  toFun q
-/

/-
Copyright (c) 2025 Anthony Wang. All rights reserved.
Released under Apache 2.0 license as described in the file LICENSE.
Authors: Anthony Wang
-/
import Mathlib


/-!
# Spherical coordinates

We define spherical coordinates similarly to polar coordinates, as an open partial homeomorphism
in `ℝ^3` between `ℝ^3 - (-∞, 0]` and `(0, +∞) × (-π, π) × (0, π)`. Its inverse is given by
`(r, θ, φ) ↦ (r sin φ cos θ, r sin φ sin θ, r cos φ)`. TODO

It satisfies the following change of variables formula (see `integral_comp_sphericalCoord_symm`):
`∫ p in sphericalCoord.target, p.1 • f (sphericalCoord.symm p) = ∫ p, f p`

-/

theorem arg_of_im_pos' {z : ℂ} (hz : 0 < z.im) : 0 < Complex.arg z := by
  have h₀ : z ≠ 0 := Ne.symm (mt (congr_arg Complex.im) hz.ne)
  suffices 0 < Real.sin z.arg by
    contrapose! this
    exact Real.sin_nonpos_of_nonnpos_of_neg_pi_le this (by linarith [Complex.neg_pi_lt_arg z])
  simp only [Complex.sin_arg, lt_div_iff₀ (norm_pos_iff.mpr h₀), zero_mul, hz]

noncomputable section Real

open Real Set MeasureTheory

open scoped ENNReal Real Topology

/-- The spherical coordinates are an open partial homeomorphism in `ℝ^3`, mapping `(r sin φ cos θ, r sin φ sin θ, r cos φ)` to
`(r, θ, φ)`. It is a homeomorphism between `ℝ^3 - (-∞, 0]` and `(0, +∞) × (-π, π) × (0, π)`. TODO -/
def sphericalCoord : PartialHomeomorph (ℝ × ℝ × ℝ) (ℝ × ℝ × ℝ) where
  toFun q := (√(q.1 ^ 2 + q.2.1 ^ 2 + q.2.2 ^ 2), Complex.arg (Complex.equivRealProd.symm (q.1, q.2.1)), Complex.arg (Complex.equivRealProd.symm (q.2.2, √(q.1 ^ 2 + q.2.1 ^ 2))))
  invFun p := (p.1 * sin p.2.2 * cos p.2.1, p.1 * sin p.2.2 * sin p.2.1, p.1 * cos p.2.2)
  source := {q | 0 < q.1} ∪ {q | q.2.1 ≠ 0}
  target := Ioi (0 : ℝ) ×ˢ Ioo (-π) π ×ˢ Ioo 0 π
  map_target' := by
    rintro ⟨r, θ, φ⟩ ⟨hr, hθ, hφ⟩
    dsimp at hr hθ hφ
    rcases eq_or_ne θ 0 with (rfl | h'θ)
    · simp only [ne_eq, cos_zero, mul_one, sin_zero, mul_zero, mem_union, mem_setOf_eq]
      left
      exact Left.mul_pos hr (sin_pos_of_mem_Ioo hφ)
    · simp only [ne_eq, mem_union, mem_setOf_eq, mul_eq_zero, not_or]
      right
      and_intros
      · linarith [mem_Ioi.mp hr]
      · linarith [sin_pos_of_mem_Ioo hφ]
      · simp [sin_eq_zero_iff_of_lt_of_lt hθ.1 hθ.2, h'θ]
  map_source' := by
    rintro ⟨x, y, z⟩ hxy
    simp only [prodMk_mem_set_prod_eq, mem_Ioi, sqrt_pos, mem_Ioo, Complex.neg_pi_lt_arg,
      true_and, Complex.arg_lt_pi_iff]
    simp
    have hpos : 0 < x ^ 2 + y ^ 2 := by
      rcases hxy with hxy | hxy
      · dsimp at hxy; linarith [sq_pos_of_ne_zero hxy.ne', sq_nonneg y]
      · linarith [sq_nonneg x, sq_pos_of_ne_zero hxy]
    and_intros
    · rcases hxy with hxy | hxy
      · dsimp at hxy; linarith [sq_pos_of_ne_zero hxy.ne', sq_nonneg y, sq_nonneg z]
      · linarith [sq_nonneg x, sq_pos_of_ne_zero hxy, sq_nonneg z]
    · rcases hxy with hxy | hxy
      · exact Or.inl (le_of_lt hxy)
      · exact Or.inr hxy
    · apply arg_of_im_pos'
      simp [hpos]
    · exact Or.inr (by simp only [sqrt_ne_zero', hpos])
  right_inv' := by
    rintro ⟨r, θ, φ⟩ ⟨hr, hθ, hφ⟩
    ext <;> dsimp at hr hθ hφ ⊢
    · conv_rhs => rw [← sqrt_sq (le_of_lt hr), ← one_mul (r ^ 2), ← sin_sq_add_cos_sq φ, ← one_mul (sin φ ^ 2), ← sin_sq_add_cos_sq θ]
      congr 1
      ring
    · convert Complex.arg_mul_cos_add_sin_mul_I (Left.mul_pos hr (sin_pos_of_mem_Ioo hφ)) ⟨hθ.1, hθ.2.le⟩
      simp only [Complex.equivRealProd_symm_apply, Complex.ofReal_mul, Complex.ofReal_sin,
        Complex.ofReal_cos]
      ring
    · have : -π < φ := by linarith [pi_pos, hφ.1]
      convert Complex.arg_mul_cos_add_sin_mul_I hr ⟨this, hφ.2.le⟩
      have : √((r * sin φ * cos θ) ^ 2 + (r * sin φ * sin θ) ^ 2) = r * sin φ := by
        conv_rhs => rw [← sqrt_sq (le_of_lt (Left.mul_pos hr (sin_pos_of_mem_Ioo hφ))), ← one_mul ((r * sin φ) ^ 2), ← sin_sq_add_cos_sq θ]
        congr 1
        ring
      simp only [this, Complex.equivRealProd_symm_apply, Complex.ofReal_mul, Complex.ofReal_cos,
        Complex.ofReal_sin]
      ring
  left_inv' := by
    rintro ⟨x, y, z⟩ _
    -- have A : √(x ^ 2 + y ^ 2) = ‖x + y * Complex.I‖ := by
    --   rw [Complex.norm_def, Complex.normSq_add_mul_I]
    -- have B : √(x ^ 2 + y ^ 2 + z ^ 2) = ‖z + √(x ^ 2 + y ^ 2) * Complex.I‖ := by
    --   rw [Complex.norm_def, Complex.normSq_add_mul_I]
    simp -- [Complex.sin_arg]
    and_intros
    ·
      rw [ Complex.sin_arg, Complex.cos_arg ] <;> aesop;
      · -- Simplify the expression by canceling out the common terms.
        field_simp [Complex.normSq, Complex.norm_def]
        ring;
        field_simp [mul_comm, mul_assoc, mul_left_comm];
      · norm_num [ Complex.normSq, Complex.norm_def ] ; ring;
        -- By simplifying, we can see that the expression indeed equals $x$.
        field_simp
        ring;
      · -- From the equation $x + y * Complex.I = 0$, we can separate the real and imaginary parts to get $x = 0$ and $y = 0$.
        have h_real_imag : x = 0 ∧ y = 0 := by
          -- Since $x + y * Complex.I = 0$, we can separate the real and imaginary parts to get $x = 0$ and $y = 0$.
          simp [Complex.ext_iff] at a_1; exact ⟨by linarith, by linarith⟩;
        -- Since $x = 0$ and $0 < x$, this is a contradiction.
        aesop;
      · simp_all +decide [ Complex.ext_iff ]
    ·
      rw [ Complex.sin_arg, Complex.sin_arg ] <;> aesop;
      · -- Simplify the norms and cancel out the common terms.
        field_simp [Complex.normSq, Complex.norm_def]
        ring;
        -- By simplifying, we can see that the expression indeed equals y.
        field_simp [mul_comm, mul_assoc, mul_left_comm];
      · -- Simplify the norms and square roots in the expression.
        field_simp [Complex.normSq, Complex.norm_def]
        ring;
        simp +decide [ mul_assoc, mul_comm, mul_left_comm, ne_of_gt ( Real.sqrt_pos.mpr ( show 0 < x ^ 2 + y ^ 2 + z ^ 2 by positivity ) ), ne_of_gt ( Real.sqrt_pos.mpr ( show 0 < x ^ 2 + y ^ 2 by positivity ) ) ]
    ·
      rw [ Complex.cos_arg ] <;> aesop;
      · norm_num [ Complex.normSq, Complex.norm_def ];
        rw [ Real.mul_self_sqrt ( by positivity ) ] ; ring_nf;
        rw [ mul_assoc, mul_inv_cancel₀ ( by positivity ), mul_one ];
      · norm_num [ Complex.normSq, Complex.norm_def ];
        rw [ Real.mul_self_sqrt ( by positivity ) ] ; ring;
        rw [ mul_assoc, mul_inv_cancel₀ ( by positivity ), mul_one ];
      · -- Since $z$ is real and $\sqrt{x^2 + y^2}$ is real, the only way their sum can be zero is if both parts are zero. However, $\sqrt{x^2 + y^2}$ is always non-negative and since $x$ is positive, it can't be zero. Therefore, $z$ must be zero, but then the imaginary part would be $\sqrt{x^2 + y^2}$, which can't be zero either. This leads to a contradiction.
        have h_contra : z = 0 ∧ Real.sqrt (x^2 + y^2) = 0 := by
          simp_all +decide [ Complex.ext_iff ];
        exact absurd h_contra.2 ( by positivity );
      · simp_all +decide [ Complex.ext_iff ];
        exact ne_of_gt ( Real.sqrt_pos.mpr ( by positivity ) ) a_1.2
  open_target := isOpen_Ioi.prod (isOpen_Ioo.prod isOpen_Ioo)
  open_source :=
    (isOpen_lt continuous_const continuous_fst).union
      (isOpen_ne_fun (Continuous.fst continuous_snd) continuous_const)
  continuousOn_invFun := by fun_prop
  continuousOn_toFun := by
    refine' ContinuousOn.prod _ _;
    · exact Continuous.continuousOn ( Real.continuous_sqrt.comp <| Continuous.add ( Continuous.add ( continuous_fst.pow 2 ) ( continuous_snd.fst.pow 2 ) ) ( continuous_snd.snd.pow 2 ) );
    · field_simp;
      -- The argument function is continuous except at zero, and since the inputs are never zero in the given set, the argument function is continuous.
      have h_arg_cont : ContinuousOn (fun q : ℝ × ℝ × ℝ => Complex.arg (q.1 + q.2.1 * Complex.I)) ({q : ℝ × ℝ × ℝ | 0 < q.1} ∪ {q : ℝ × ℝ × ℝ | ¬q.2.1 = 0}) ∧ ContinuousOn (fun q : ℝ × ℝ × ℝ => Complex.arg (q.2.2 + Real.sqrt (q.1^2 + q.2.1^2) * Complex.I)) ({q : ℝ × ℝ × ℝ | 0 < q.1} ∪ {q : ℝ × ℝ × ℝ | ¬q.2.1 = 0}) := by
        constructor <;> refine' ContinuousAt.continuousOn _ <;> intro q hq <;> aesop;
        · refine' Complex.continuousAt_arg _ |> fun h => h.comp _;
          · simp [Complex.slitPlane, h];
          · fun_prop (disch := norm_num);
        · refine' Complex.continuousAt_arg _ |> fun h => h.comp _;
          · simp [Complex.slitPlane, h_1];
          · fun_prop (disch := norm_num);
        · refine' Complex.continuousAt_arg _ |> ContinuousAt.comp <| _;
          · norm_num [ Complex.slitPlane, Complex.ext_iff ];
            exact Or.inr ( ne_of_gt ( Real.sqrt_pos.mpr ( by positivity ) ) );
          · fun_prop;
        · refine' Complex.continuousAt_arg _ |> fun h => h.comp _;
          · norm_num [ Complex.slitPlane ];
            exact Or.inr ( ne_of_gt ( Real.sqrt_pos.mpr ( by positivity ) ) );
          · fun_prop (disch := norm_num);
      exact ContinuousOn.prod h_arg_cont.1 h_arg_cont.2