Mercurial > repos > public > sbplib_julia

comparison src/Grids/parameter_space.jl @ 2057:8a2a0d678d6f feature/lazy_tensors/pretty_printing

Find changesets by keywords (author, files, the commit message), revision number or hash, or revset expression.
Merge default
author Jonatan Werpers <jonatan@werpers.com>
date Tue, 10 Feb 2026 22:41:19 +0100
parents 889c18ad56bf
children
comparison
equal deleted inserted replaced
1110:c0bff9f6e0fb 2057:8a2a0d678d6f
1 """
2 ParameterSpace{D}
3
4 A space of parameters of dimension `D`.
5
6 Common parameter spaces are created using functions for unit sized spaces
7 * [`unitinterval`](@ref)
8 * [`unitsquare`](@ref)
9 * [`unitcube`](@ref)
10 * [`unithyperbox`](@ref)
11 * [`unittriangle`](@ref)
12 * [`unittetrahedron`](@ref)
13 * [`unitsimplex`](@ref)
14
15 See also: [`Interval`](@ref), [`HyperBox`](@ref),
16 [`Simplex`](@ref).
17 """
18 abstract type ParameterSpace{D} end
19 Base.ndims(::ParameterSpace{D}) where D = D
20
21 @doc """
22 in(x, S::ParameterSpace)
23 ∈(x, S::ParameterSpace)
24
25 Test if the point `x` is in the parameter space `S`.
26 """ Base.in(x,::ParameterSpace)
27
28 """
29 Interval{T} <: ParameterSpace{1}
30
31 A `ParameterSpace` representing an interval.
32 """
33 struct Interval{T} <: ParameterSpace{1}
34 a::T
35 b::T
36 end
37
38 """
39 Interval(a,b)
40
41 An interval with limits `a` and `b`.
42 """
43 function Interval(a,b)
44 a, b = promote(a, b)
45 Interval{typeof(a)}(a,b)
46 end
47
48 """
49 limits(i::Interval)
50
51 The limits of the interval.
52 """
53 limits(i::Interval) = (i.a, i.b)
54
55 boundary_identifiers(::Interval) = (LowerBoundary(), UpperBoundary())
56
57 Base.in(x, i::Interval) = i.a <= x <= i.b
58
59 """
60 unitinterval(T=Float64)
61
62 The interval ``(0,1)``.
63 """
64 unitinterval(T=Float64) = Interval(zero(T), one(T))
65
66
67 """
68 HyperBox{T,D} <: ParameterSpace{D}
69
70 A `ParameterSpace` representing a hyper box.
71 """
72 struct HyperBox{T,D} <: ParameterSpace{D}
73 a::SVector{D,T}
74 b::SVector{D,T}
75 end
76
77 """
78 HyperBox(a,b)
79
80 A `HyperBox` with lower limits `a` and upper limits `b` for each dimension.
81 """
82 function HyperBox(a,b)
83 ET = promote_type(eltype(a),eltype(b))
84 T = SVector{length(a),ET}
85 HyperBox(convert(T,a), convert(T,b))
86 end
87
88 Rectangle{T} = HyperBox{T,2}
89 Box{T} = HyperBox{T,3}
90
91 """
92 limits(box::HyperBox, d)
93
94 Limits of `box` along dimension `d`.
95 """
96 limits(box::HyperBox, d) = (box.a[d], box.b[d])
97
98 """
99 limits(box::HyperBox)
100
101 The lower and upper limits of `box` as tuples.
102 """
103 limits(box::HyperBox) = (box.a, box.b)
104
105 function boundary_identifiers(box::HyperBox)
106 mapreduce(vcat, 1:ndims(box)) do d
107 [
108 CartesianBoundary{d, LowerBoundary}(),
109 CartesianBoundary{d, UpperBoundary}(),
110 ]
111 end
112 end
113
114 function Base.in(x, box::HyperBox)
115 return all(eachindex(x)) do i
116 box.a[i] <= x[i] <= box.b[i]
117 end
118 end
119
120 """
121 unitsquare(T=Float64)
122
123 The square limited by 0 and 1 in each dimension.
124 """
125 unitsquare(T=Float64) = unithyperbox(T,2)
126
127 """
128 unitcube(T=Float64)
129
130 The cube limited by 0 and 1 in each dimension.
131 """
132 unitcube(T=Float64) = unithyperbox(T,3)
133
134 """
135 unithyperbox(T=Float64, D)
136
137 The hypercube limited by 0 and 1 in each dimension.
138 """
139 unithyperbox(T, D) = HyperBox((@SVector zeros(T,D)), (@SVector ones(T,D)))
140 unithyperbox(D) = unithyperbox(Float64,D)
141
142
143 """
144 Simplex{T,D,NV} <: ParameterSpace{D}
145
146 A `ParameterSpace` representing a simplex.
147 """
148 struct Simplex{T,D,NV} <: ParameterSpace{D}
149 verticies::NTuple{NV,SVector{D,T}}
150
151 Simplex(verticies::Tuple{SVector{D,T}, Vararg{SVector{D,T},N}}) where {T,D,N} = new{T,D,N+1}(verticies)
152 Simplex(::Tuple{}) = throw(ArgumentError("Must provide at least one vertex."))
153 end
154
155 """
156 Simplex(verticies...)
157
158 A simplex with the given verticies.
159 """
160 function Simplex(verticies::Vararg{AbstractArray})
161 ET = mapreduce(eltype,promote_type,verticies)
162 T = SVector{length(verticies[1]),ET}
163
164 return Simplex(Tuple(convert(T,v) for v ∈ verticies))
165 end
166
167 function Base.in(x, s::Simplex)
168 v₁ = s.verticies[1]
169 V = map(s.verticies) do v
170 v - v₁
171 end
172
173 A = hcat(V[2:end]...) # Matrix with edge vectors as columns
174 λ = A \ (x - v₁)
175
176 λ_full = (1 - sum(λ), λ...) # Full barycentric coordinates
177
178 return all(λᵢ -> zero(λᵢ) ≤ λᵢ ≤ one(λᵢ), λ_full)
179 end
180
181 """
182 verticies(s::Simplex)
183
184 Verticies of `s`.
185 """
186 verticies(s::Simplex) = s.verticies
187
188 Triangle{T} = Simplex{T,2}
189 Tetrahedron{T} = Simplex{T,3}
190
191 """
192 unittriangle(T=Float64)
193
194 The simplex with verticies ``(0,0)``, ``(1,0)``, and ``(0,1)``.
195 """
196 unittriangle(T=Float64) = unitsimplex(T,2)
197
198 """
199 unittetrahedron(T=Float64)
200
201 The simplex with verticies ``(0,0,0)``, ``(1,0,0)``, ``(0,1,0)``, and ``(0,0,1)``.
202 """
203 unittetrahedron(T=Float64) = unitsimplex(T,3)
204
205 """
206 unitsimplex(T=Float64,D)
207
208 The unit simplex in dimension `D` with verticies ``(0,0,0,...)``, ``(1,0,0,...)``, ``(0,1,0,...)``, ``(0,0,1,...)``...
209 """
210 function unitsimplex(T,D)
211 z = @SVector zeros(T,D)
212 unitelement = one(eltype(z))
213 verticies = ntuple(i->setindex(z, unitelement, i), D)
214 return Simplex((z,verticies...))
215 end
216 unitsimplex(D) = unitsimplex(Float64, D)