changeset 2012:4617e4b74b82 feature/grids/geometry_functions

Add arc() for constructing circle arcs between two points
author Jonatan Werpers <jonatan@werpers.com>
date Wed, 07 May 2025 08:39:06 +0200
parents d0b6c63c506e
children 7895b509f9bf
files src/Grids/geometry.jl test/Grids/geometry_test.jl
diffstat 2 files changed, 92 insertions(+), 0 deletions(-) [+]
line wrap: on
line diff
diff -r d0b6c63c506e -r 4617e4b74b82 src/Grids/geometry.jl
--- a/src/Grids/geometry.jl	Wed May 07 08:38:35 2025 +0200
+++ b/src/Grids/geometry.jl	Wed May 07 08:39:06 2025 +0200
@@ -181,6 +181,40 @@
     return nothing
 end
 
+
+"""
+    arc(a,b,r)
+
+# TODO
+"""
+function arc(a,b,r)
+    if abs(r) < norm(b-a)/2
+        throw(DomainError(r, "arc was called with radius r = $r smaller than half the distance between the points."))
+    end
+
+    R̂ = @SMatrix[0 -1; 1 0]
+
+    α = sign(r)*√(r^2 - norm((b-a)/2)^2)
+    t̂ = R̂*(b-a)/norm(b-a)
+
+    c = (a+b)/2 + α*t̂
+
+    ca = a-c
+    cb = b-c
+    θₐ = atan(ca[2],ca[1])
+    θᵦ = atan(cb[2],cb[1])
+
+    Δθ = mod(θᵦ-θₐ+π, 2π)-π # Δθ in the interval (-π,π)
+
+    if r > 0
+        Δθ = abs(Δθ)
+    else
+        Δθ = -abs(Δθ)
+    end
+
+    return Arc(Circle(c,abs(r)), θₐ, θₐ+Δθ)
+end
+
 """
     TransfiniteInterpolationSurface(c₁, c₂, c₃, c₄)
 
diff -r d0b6c63c506e -r 4617e4b74b82 test/Grids/geometry_test.jl
--- a/test/Grids/geometry_test.jl	Wed May 07 08:38:35 2025 +0200
+++ b/test/Grids/geometry_test.jl	Wed May 07 08:39:06 2025 +0200
@@ -181,6 +181,64 @@
     end
 end
 
+@testset "arc" begin
+    a = [0,0]
+    b = [1,0]
+    A = arc(a,b,1/2)
+    @test A(0) ≈ a atol=1e-15
+    @test A(1) ≈ b
+    @test A(0.5) ≈ [0.5, -0.5]
+
+    A = arc(a,b,-1/2)
+    @test A(0) ≈ a atol=1e-15
+    @test A(1) ≈ b
+    @test A(0.5) ≈ [0.5, 0.5]
+
+    @testset "Unit arc" begin
+        A = arc([1,0],[0,1],1)
+        @test A(0) ≈ [1,0]
+        @test A(1) ≈ [0,1]
+        @testset for t ∈ range(0,1,13)
+            @test A(t) ≈ [cos(t*π/2), sin(t*π/2)]
+        end
+    end
+
+    @testset "Inverted unit arc" begin
+        A = arc([1,0],[0,1],-1)
+        @test A(0) ≈ [1,0]
+        @test A(1) ≈ [0,1]
+        @testset "Inverted unit arc t=$t" for t ∈ range(0,1,13)
+            @test A(t) ≈ [1+cos(-π/2 - t*π/2), 1+sin(-π/2 - t*π/2)]
+        end
+    end
+
+    @testset "Quarters of unit circle" begin
+        unitvec(θ) = [cos(θ), sin(θ)]
+        @testset "θ ∈ ($(i)π/4, $(i+2)π/4)" for i ∈ range(0, step=1, length=8)
+            θ = i*π/4
+            @testset let θ₀ = θ, θ₁ = θ+π/2, r = 1
+                A = arc(unitvec(θ₀), unitvec(θ₁), r)
+                @test A(0) ≈ unitvec(θ)
+                @test A(1/3) ≈ unitvec(θ+π/6)
+                @test A(1/2) ≈ unitvec(θ+π/4)
+                @test A(2/3) ≈ unitvec(θ+π/3)
+                @test A(1) ≈ unitvec(θ+π/2)
+            end
+
+            @testset let θ₀ = θ+π/2, θ₁ = θ, r = -1
+                A = arc(unitvec(θ₀), unitvec(θ₁), r)
+                @test A(0) ≈ unitvec(θ+π/2)
+                @test A(1/3) ≈ unitvec(θ+π/3)
+                @test A(1/2) ≈ unitvec(θ+π/4)
+                @test A(2/3) ≈ unitvec(θ+π/6)
+                @test A(1) ≈ unitvec(θ)
+            end
+        end
+    end
+
+    @test_throws DomainError arc([-1,0], [1,0], 0.7)
+end
+
 @testset "TransfiniteInterpolationSurface" begin
     @testset "Constructors" begin
         @test TransfiniteInterpolationSurface(t->[1,2], t->[2,1], t->[0,0], t->[1,1]) isa TransfiniteInterpolationSurface