From 582ccf4f48970267e94cd531df63eeca8f816f65 Mon Sep 17 00:00:00 2001 From: Tom McGrath Date: Sun, 9 Feb 2025 15:06:04 -0600 Subject: [PATCH 1/2] Implement sin, cos `sin` and `cos` of a `GVar` return the exact first three moments of the transformed Gaussian, obtained from its characteristic function, and propagate the reliability diagnostics as `exp` does. `sincos` returns the pair. Assisted-by: Claude Opus 4.8 --- src/GaussianRandomVariables.jl | 38 +++++++++++++++++++++++++++++++++- 1 file changed, 37 insertions(+), 1 deletion(-) diff --git a/src/GaussianRandomVariables.jl b/src/GaussianRandomVariables.jl index 4c581e4..bd1712b 100644 --- a/src/GaussianRandomVariables.jl +++ b/src/GaussianRandomVariables.jl @@ -5,7 +5,7 @@ using ThickNumbers import Base: +, -, *, /, //, ^, inv import Base: abs, abs2, max, min, sqrt -import Base: log, exp +import Base: log, exp, sin, cos, sincos export GVar, ± export skewness, moment_error, distrust @@ -410,4 +410,40 @@ sqrt(a::GVar{<:AbstractFloat}) = isempty(a) ? a : gmap(sqrt, x -> 1/(2*sqrt(x)), x -> -1/(4*sqrt(x^3)), x -> 3/(8*sqrt(x^5)), x -> -15/(16*sqrt(x^7)), a) +# sin and cos of a Gaussian are not Gaussian, but their first three moments are +# exact in closed form from the characteristic function E[exp(itx)] = exp(itc - t²σ²/2). +# With s = sin c, k = cos c and the Gaussian damping factors d1 = exp(-σ²/2), +# d2 = exp(-2σ²), d9 = exp(-9σ²/2): +# +# E[sin x] = s·d1 E[cos x] = k·d1 +# E[sin²x] = (1 - cos2c·d2)/2 E[cos²x] = (1 + cos2c·d2)/2 +# E[sin³x] = (3s·d1 - sin3c·d9)/4 E[cos³x] = (3k·d1 + cos3c·d9)/4 +# +# The mean is exact in σ, so nothing is charged to `err`; as for `exp`, only the +# input's own skew and center error propagate, to leading order through the local +# derivatives (f' = k, f'' = -s, f''' = -k for sin; negated/rotated for cos). +function sin(a::GVar{<:AbstractFloat}) + isempty(a) && return a + c, σ, κ3, err = a.center, a.σ, a.κ3, a.err + d1, d2, d9 = exp(-σ^2/2), exp(-2σ^2), exp(-9σ^2/2) + s, k = sin(c), cos(c) + μ = s*d1 + v = (1 - cos(2c)*d2)/2 - μ^2 + e3 = (3s*d1 - sin(3c)*d9)/4 + return assemble(μ - k*κ3/6, v - s*k*κ3, (e3 - 3μ*(v + μ^2) + 2μ^3) + k^3*κ3, abs(k*d1)*err) +end + +function cos(a::GVar{<:AbstractFloat}) + isempty(a) && return a + c, σ, κ3, err = a.center, a.σ, a.κ3, a.err + d1, d2, d9 = exp(-σ^2/2), exp(-2σ^2), exp(-9σ^2/2) + s, k = sin(c), cos(c) + μ = k*d1 + v = (1 + cos(2c)*d2)/2 - μ^2 + e3 = (3k*d1 + cos(3c)*d9)/4 + return assemble(μ + s*κ3/6, v + s*k*κ3, (e3 - 3μ*(v + μ^2) + 2μ^3) - s^3*κ3, abs(s*d1)*err) +end + +sincos(a::GVar) = (sin(a), cos(a)) + end # module From 6badbf78a97dce17adf14d163a741a189ebcb749 Mon Sep 17 00:00:00 2001 From: Tim Holy Date: Wed, 15 Jul 2026 00:38:37 -0500 Subject: [PATCH 2/2] Test sin and cos Assisted-by: Claude Opus 4.8 --- test/runtests.jl | 46 ++++++++++++++++++++++++++++++++++++++++++++++ 1 file changed, 46 insertions(+) diff --git a/test/runtests.jl b/test/runtests.jl index d962183..8a5c4b8 100644 --- a/test/runtests.jl +++ b/test/runtests.jl @@ -193,6 +193,52 @@ ispositive(x) = x > 0 @test abs2(3 ± 0.5) ⩪ (3 ± 0.5)^2 end + # sin and cos of a Gaussian are not Gaussian, but their first three moments are + # exact, so mid, rad and κ3 match closed forms and nothing is charged to `err`. + @testset "sin and cos" begin + c, σ = 0.7, 0.5 + d1, d2 = exp(-σ^2/2), exp(-2σ^2) + s = sin(GVar(c, σ)) + @test mid(s) ≈ sin(c)*d1 + @test rad(s)^2 ≈ (1 - cos(2c)*d2)/2 - (sin(c)*d1)^2 + @test moment_error(s) == 0 + k = cos(GVar(c, σ)) + @test mid(k) ≈ cos(c)*d1 + @test rad(k)^2 ≈ (1 + cos(2c)*d2)/2 - (cos(c)*d1)^2 + @test moment_error(k) == 0 + + # A zero-mean input is symmetric: sin has zero mean and no skew, while cos + # peaks at 1 and folds its mass downward, so it is skewed toward low values. + @test mid(sin(0 ± σ)) == 0 + @test sin(0 ± σ).κ3 == 0 + @test mid(cos(0 ± σ)) ≈ exp(-σ^2/2) + @test skewness(cos(0 ± σ)) < 0 + + # A vanishing spread is an ordinary number. + @test sin(GVar(c, 0.0)) ⩪ GVar(sin(c), 0.0) + @test cos(GVar(c, 0.0)) ⩪ GVar(cos(c), 0.0) + + # `sincos` returns the pair; the empty set propagates through all three. + sc = sincos(GVar(c, σ)) + @test sc[1] ⩪ sin(GVar(c, σ)) && sc[2] ⩪ cos(GVar(c, σ)) + e = GVar(0.0, -1.0) + @test isempty(sin(e)) && isempty(cos(e)) + @test all(isempty, sincos(e)) + + @test @inferred(sin(GVar(1.0, 0.5))) isa GVar{Float64} + @test @inferred(sincos(GVar(1.0, 0.5))) isa Tuple{GVar{Float64},GVar{Float64}} + + # Sampled moments confirm the closed forms across a range of spreads. + for (μ, σ) in ((0.6, 0.2), (1.0, 0.5), (2.0, 0.8), (0.9, 1.0)) + @test testscalar(sin, μ, σ; rtol=0.02, n=10^6) + @test testscalar(cos, μ, σ; rtol=0.02, n=10^6) + end + for (μ, σ) in ((1.0, 0.5), (2.0, 0.8), (0.9, 1.0)) + @test testskew(sin, μ, σ) + @test testskew(cos, μ, σ) + end + end + @testset "min and max" begin a, b = GVar(3.0, 1.0), GVar(3.5, 1.0) # spans [2,4] and [2.5,4.5] @test loval(min(a, b)) ≈ 2.0 && hival(min(a, b)) ≈ 4.0