-- Quantitative descent along a complex direction, without root extraction.
-- Part of Cumulus, the Grace17 mathematics library.
-- Copyright (c) 2026 Grant Passmore.
-- Released under the Apache License 2.0. See LICENSE in the repository root.
-- Authors: Grant Passmore

module Grace.Cumulus.Complex.Descent
import LocalPolynomial
import Real.Power

def complexScale (t : Real) (z : Euclidean.Point) : Set :=
  (t * Euclidean.x(z), t * Euclidean.y(z))

section Complexes
notation ring ℂ := (complexCarrier, complexAdd, complexZero, complexNeg, complexMul, complexOne)
notation ring R := (ℝ, realAdd, realZero, realNeg, realMul, realOne)

theorem complexScale.closed (t : Real) (z : Set)
  requires member : z ∈ ℂ
  shows complexScale(t, z) ∈ ℂ
proof
  have point := complex.point with *
  unfold complexScale
  apply complexCarrier.pair

theorem complexScale.mul (t : Real) (z : Set)
  requires member : z ∈ ℂ
  shows complexScale(t, z) = ℂ[(t, realZero) * z]
proof
  obtain ⟨a, aReal, b, bReal, shape⟩ := complexCarrier.cases at z with *
  rw shape
  unfold complexScale
  simp [Euclidean.x.pair, Euclidean.y.pair, complexMul.pair]
  have firstCoordinate : R[t * a - 0 * b] = R[t * a] by ring
  have secondCoordinate : R[t * b + 0 * a] = R[t * b] by ring
  rw [firstCoordinate, secondCoordinate]
  rfl

theorem complexNormSquared.scale (t : Real) (z : Set)
  requires member : z ∈ ℂ
  shows complexNormSquared(complexScale(t, z)) = R[t * t * complexNormSquared(z)]
proof
  obtain ⟨a, aReal, b, bReal, shape⟩ := complexCarrier.cases at z with *
  rw shape
  unfold complexScale
  simp [Euclidean.x.pair, Euclidean.y.pair, complexNormSquared.pair]
  ring

theorem complexScale.dot_right (t : Real) (z w : Set)
  requires first : z ∈ ℂ
  requires second : w ∈ ℂ
  shows Euclidean.dot(z, complexScale(t, w)) = R[t * Euclidean.dot(z, w)]
proof
  obtain ⟨a, aReal, b, bReal, zShape⟩ := complexCarrier.cases at z with *
  obtain ⟨c, cReal, d, dReal, wShape⟩ := complexCarrier.cases at w with *
  rw [zShape, wShape]
  unfold complexScale
  simp [Euclidean.dot, Euclidean.x.pair, Euclidean.y.pair]
  ring

theorem complexScale.product (t s : Real) (z w : Set)
  requires first : z ∈ ℂ
  requires second : w ∈ ℂ
  shows ℂ[complexScale(t, z) * complexScale(s, w)] = complexScale(R[t * s], ℂ[z * w])
proof
  obtain ⟨a, aReal, b, bReal, zShape⟩ := complexCarrier.cases at z with *
  obtain ⟨c, cReal, d, dReal, wShape⟩ := complexCarrier.cases at w with *
  rw [zShape, wShape]
  unfold complexScale
  simp [complexMul.pair, Euclidean.x.pair, Euclidean.y.pair]
  have firstCoordinate : R[(t * a) * (s * c) - (t * b) * (s * d)] = R[(t * s) * (a * c - b * d)] by ring
  have secondCoordinate : R[(t * a) * (s * d) + (t * b) * (s * c)] = R[(t * s) * (a * d + b * c)] by ring
  rw [firstCoordinate, secondCoordinate]
  rfl

theorem complexScale.one (z : Set)
  requires member : z ∈ ℂ
  shows complexScale(realOne, z) = z
proof
  obtain ⟨a, aReal, b, bReal, shape⟩ := complexCarrier.cases at z with *
  rw shape
  unfold complexScale
  simp [Euclidean.x.pair, Euclidean.y.pair, realMulValue.one_left]

theorem complexScale.mul_right (t : Real) (z w : Set)
  requires first : z ∈ ℂ
  requires second : w ∈ ℂ
  shows ℂ[complexScale(t, z) * w] = complexScale(t, ℂ[z * w])
proof
  have product := complexScale.product at t, realOne, z, w with *
  rwa [complexScale.one, realMulValue.one_right] at product

theorem complexScale.power (t : Real) (z : Set) (n : Nat)
  requires member : z ∈ ℂ
  shows ℂ[complexScale(t, z) ^ n] = complexScale(realPow(t, n), ℂ[z ^ n])
proof
  induction n
  | zero =>
    rw [monoidPow.zero, monoidPow.zero, realPow.zero]
    unfold complexScale
    unfold complexOne
    simp [Euclidean.x.pair, Euclidean.y.pair, realMulValue.one_left]
  | succ(k) using ih =>
    have ring := complexes.comm_ring
    have ordinary := IsCommRing.ring with ring
    have monoid := IsRing.multiplicative with ordinary
    rw [monoidPow.succ, monoidPow.succ, realPow.succ, ih]
    exact complexScale.product at t, realPow(t, k), z, ℂ[z ^ k] with *

-- The threshold is explicit: half(-dot(z,w))/(abs(normSquared(w))+1).
theorem complexNormSquared.directional_descent (z w : Set) (t : Real)
  requires first : z ∈ ℂ
  requires second : w ∈ ℂ
  requires direction : Euclidean.dot(z, w) < realZero
  requires positive : realZero < t
  requires small : t < scaledRadius(R[-Euclidean.dot(z, w)], limitScale(complexNormSquared(w)))
  shows complexNormSquared(ℂ[z + complexScale(t, w)]) < complexNormSquared(z)
proof
  have pointZ := complex.point at z with *
  have pointW := complex.point at w with *
  have epsilon := realInv.negative_argument_positive at Euclidean.dot(z, w) with *
  have nonnegative : realZero ≤ t
    exact realLess.implies_le with *
  have control := scaledRadius.product_control with epsilon, nonnegative
  have shifted := realLess.add_left at R[t * complexNormSquared(w)], R[-Euclidean.dot(z, w)], Euclidean.dot(z, w) with *
  have cancel : R[Euclidean.dot(z, w) + -Euclidean.dot(z, w)] = realZero by ring
  rw cancel at shifted
  have combined := realLess.add at Euclidean.dot(z, w), realZero, R[Euclidean.dot(z, w) + t * complexNormSquared(w)], realZero with *
  rw realAddValue.zero_left at combined
  have normalize : R[Euclidean.dot(z, w) + (Euclidean.dot(z, w) + t * complexNormSquared(w))] = R[2 * Euclidean.dot(z, w) + t * complexNormSquared(w)]
    ring
  rw normalize at combined
  have negative : R[2 * Euclidean.dot(z, w) + t * complexNormSquared(w)] < realZero by exact combined
  have multiplied := realLess.mul_left_of_positive at R[2 * Euclidean.dot(z, w) + t * complexNormSquared(w)], realZero, t with *
  rw realMulValue.zero_right at multiplied
  have descent : R[2 * Euclidean.dot(z, complexScale(t, w)) + complexNormSquared(complexScale(t, w))] < realZero
    rw [complexScale.dot_right, complexNormSquared.scale]
    have regroup : R[2 * (t * Euclidean.dot(z, w)) + t * t * complexNormSquared(w)] = R[t * (2 * Euclidean.dot(z, w) + t * complexNormSquared(w))]
      ring
    rwa regroup
  exact complexNormSquared.strict_descent with *

theorem complexNormSquared.descent_radius (z w : Set)
  requires first : z ∈ ℂ
  requires second : w ∈ ℂ
  requires direction : Euclidean.dot(z, w) < realZero
  shows ∃ d ∈ ℝ, realZero < d ∧
    (∀ t ∈ ℝ, realZero < t → t < d → complexNormSquared(ℂ[z + complexScale(t, w)]) < complexNormSquared(z))
proof
  have pointZ := complex.point at z with *
  have pointW := complex.point at w with *
  have epsilon := realInv.negative_argument_positive with *
  have scalePositive := limitScale.positive at complexNormSquared(w) with *
  have radiusIn : scaledRadius(R[-Euclidean.dot(z, w)], limitScale(complexNormSquared(w))) ∈ ℝ by typecheck
  have radiusPositive := scaledRadius.positive at R[-Euclidean.dot(z, w)], limitScale(complexNormSquared(w)) with *
  refine ⟨scaledRadius(R[-Euclidean.dot(z, w)], limitScale(complexNormSquared(w))), radiusIn, radiusPositive, ?_⟩
  intros t tIn positive small
  exact complexNormSquared.directional_descent with *

theorem complexNormSquared.descent_step (z w : Set)
  requires first : z ∈ ℂ
  requires second : w ∈ ℂ
  requires direction : Euclidean.dot(z, w) < realZero
  shows ∃ t ∈ ℝ, realZero < t ∧ complexNormSquared(ℂ[z + complexScale(t, w)]) < complexNormSquared(z)
proof
  obtain ⟨d, dIn, positive, decreases⟩ := complexNormSquared.descent_radius at z, w with *
  have halfIn : half(d) ∈ ℝ by typecheck
  have halfPositive := half.positive with *
  have halfLess := half.lt_self with positive
  have smaller := decreases with half(d), halfPositive, halfLess
  refine ⟨half(d), halfIn, halfPositive, smaller⟩

-- The local normal form restricted to an arbitrary real-scaled direction.
theorem complexPolynomial.ray_normal_form (c a : Set) (cs : List(ℂ))
  requires coefficient : c ∈ ℂ
  requires center : a ∈ ℂ
  shows ∃ n ∈ ℕ, ∃ ds, List(ℂ, ds) ∧ complexPolynomial(ds)(ℂ[0]) ≠ ℂ[0] ∧
    length{ℂ}(cs) = Grace.Core.Nat.Rec.add(n, length{ℂ}(ds)) ∧
    ComplexContinuous(complexPolynomial(ds)) ∧
    (∀ u ∈ ℂ, ∀ t ∈ ℝ, complexPolynomial(c :: cs)(ℂ[a + complexScale(t, u)]) = ℂ[complexPolynomial(c :: cs)(a) + complexScale(realPow(t, succ(n)), ℂ[u ^ succ(n) * complexPolynomial(ds)(complexScale(t, u))])])
proof
  obtain ⟨n, nIn, ds, dsList, nonzero, degree, continuous, expansion⟩ := complexPolynomial.local_normal_form at c, a, cs with *
  refine ⟨n, nIn, ds, dsList, nonzero, degree, continuous, ?_⟩
  intros u uIn t tIn
  have nNat := nat.omega.backward at n with *
  have ring := complexes.comm_ring
  have ordinary := IsCommRing.ring with ring
  have monoid := IsRing.multiplicative with ordinary
  have scaled := complexScale.closed at t, u with *
  have power := monoidPow.closed at ℂ, complexMul, complexOne, u, succ(n) with *
  have function := complexPolynomial.function at ds with *
  have value := FunctionBetween.apply_closed at complexPolynomial(ds), ℂ, ℂ, complexScale(t, u) with *
  have equation := expansion with complexScale(t, u)
  rw [equation, complexScale.power, complexScale.mul_right]
  rfl

end Complexes