Reflector antenna problem Boris Thibert LJK Universit e de - - PowerPoint PPT Presentation
Reflector antenna problem Boris Thibert LJK Universit e de Grenoble Joint work with Quentin M erigot and Pedro Machado Modelisation with optimal transport October 3-4, 2013 1 Far-Field Reflector Antenna Problem Punctual light at origin
1
Modelisation with optimal transport October 3-4, 2013
LJK Universit´ e de Grenoble Joint work with Quentin M´ erigot and Pedro Machado
2
S2
1
Punctual light at origin o, µ measure on S2
1
S2
⌫
2
S2
1
Punctual light at origin o, µ measure on S2
1
S2
(S1, ⌫) under reflection by Snell’s law.
µ ⌫
2
S2
1
Punctual light at origin o, µ measure on S2
1
S2
(S1, ⌫) under reflection by Snell’s law.
I R is parameterized over S2
⌫
2
S2
1
Punctual light at origin o, µ measure on S2
1
S2
(S1, ⌫) under reflection by Snell’s law.
I R is parameterized over S2
TR : x 2 S2
0 7! y = x 2hx|nin
µ ⌫
2
S2
1
Punctual light at origin o, µ measure on S2
1
S2
(S1, ⌫) under reflection by Snell’s law.
I R is parameterized over S2
TR : x 2 S2
0 7! y = x 2hx|nin
Brenier formulation T]µ = ⌫ µ ⌫
2
S2
1
Punctual light at origin o, µ measure on S2
1
S2
(S1, ⌫) under reflection by Snell’s law.
I R is parameterized over S2
TR : x 2 S2
0 7! y = x 2hx|nin
Brenier formulation i.e. for every borelian B T]µ = ⌫ µ(T 1(B)) = ⌫(B) µ ⌫
2
S2
1
Punctual light at origin o, µ measure on S2
1
S2
(S1, ⌫) under reflection by Snell’s law.
I R is parameterized over S2
TR : x 2 S2
0 7! y = x 2hx|nin
Brenier formulation i.e. for every borelian B Monge-Ampere equation g(T(x)) det(DT(x)) = f(x) T]µ = ⌫ µ(T 1(B)) = ⌫(B) If µ(x) = f(x)dx and ⌫(y) = g(y)dy I highly non linear µ ⌫
2
S2
1
Punctual light at origin o, µ measure on S2
1
S2
(S1, ⌫) under reflection by Snell’s law.
I R is parameterized over S2
TR : x 2 S2
0 7! y = x 2hx|nin
Brenier formulation i.e. for every borelian B Monge-Ampere equation g(T(x)) det(DT(x)) = f(x)
Caffarelli & Oliker 94
T]µ = ⌫ µ(T 1(B)) = ⌫(B) If µ(x) = f(x)dx and ⌫(y) = g(y)dy I highly non linear I Existence
Wang 96, Guan & Wang 98
I Regularity, uniqueness µ ⌫
3
S2
1
Punctual light at origin o, µ measure on S2
3
S2
1
Punctual light at origin o, µ measure on S2
1
S2
3
S2
1
Punctual light at origin o, µ measure on S2
1
S2
y1 R : paraboloid of direction y1 and focal O
3
S2
1
y1 y2 y3 Punctual light at origin o, µ measure on S2
i ⌫iyi on S2 1
3
y1 y2 y3 Punctual light at origin o, µ measure on S2
direction yi and focal distance i R(~ ) = @
i=1Pi(i)
P2 µ P1 Prescribed far-field: ⌫ = P
i ⌫iyi on S2 1
3
y1 y2 y3 Punctual light at origin o, µ measure on S2
direction yi and focal distance i R(~ ) = @
i=1Pi(i)
) = ⇡S2
) \ @Pi(i))
R(~ ) \ @P3(3) PI3(~ )
P2 µ = directions that are reflected towards yi. Prescribed far-field: ⌫ = P
i ⌫iyi on S2 1
3
y1 y2 y3 Punctual light at origin o, µ measure on S2
direction yi and focal distance i R(~ ) = @
i=1Pi(i)
) = ⇡S2
) \ @Pi(i)) Problem (FF): Find 1, . . . , N such that for every i, µ(PIi(~ )) = ⌫i.
R(~ ) \ @P3(3) PI3(~ ) amount of light reflected in direction yi.
P2 µ = directions that are reflected towards yi. Prescribed far-field: ⌫ = P
i ⌫iyi on S2 1
4
Lemma: With c(x, y) = log(1 hx|yi), and i := log(i), PIi(~ ) = {x 2 S2
0, c(x, yi) + i c(x, yj) + j
8j}. P1
PI3(~ )
P2
Caffarelli-Oliker ’94
4
Proof: @Pi(i) is parameterized in radial coordinates by Lemma: With c(x, y) = log(1 hx|yi), and i := log(i), PIi(~ ) = {x 2 S2
0, c(x, yi) + i c(x, yj) + j
8j}. P1
PI3(~ )
P2 ⇢i : x 2 S2
i 1hx|yii
Caffarelli-Oliker ’94
4
Proof: @Pi(i) is parameterized in radial coordinates by ( ) i 2 arg minj
j 1hx|yji
x 2 PIi(~ ) Lemma: With c(x, y) = log(1 hx|yi), and i := log(i), PIi(~ ) = {x 2 S2
0, c(x, yi) + i c(x, yj) + j
8j}. P1
PI3(~ )
P2 ⇢i : x 2 S2
i 1hx|yii
Caffarelli-Oliker ’94
4
Proof: @Pi(i) is parameterized in radial coordinates by ( ) i 2 arg minj
j 1hx|yji
x 2 PIi(~ ) ( ) i 2 arg minj log(j) log(1 hx|yji) Lemma: With c(x, y) = log(1 hx|yi), and i := log(i), PIi(~ ) = {x 2 S2
0, c(x, yi) + i c(x, yj) + j
8j}. P1
PI3(~ )
P2 ⇢i : x 2 S2
i 1hx|yii
Caffarelli-Oliker ’94
4
Proof: @Pi(i) is parameterized in radial coordinates by ( ) i 2 arg minj
j 1hx|yji
x 2 PIi(~ ) ( ) i 2 arg minj log(j) log(1 hx|yji) ( ) i 2 arg minj j + c(x, yj) Lemma: With c(x, y) = log(1 hx|yi), and i := log(i), PIi(~ ) = {x 2 S2
0, c(x, yi) + i c(x, yj) + j
8j}. P1
PI3(~ )
P2 ⇢i : x 2 S2
i 1hx|yii
Caffarelli-Oliker ’94
4
Proof: @Pi(i) is parameterized in radial coordinates by ( ) i 2 arg minj
j 1hx|yji
x 2 PIi(~ ) ( ) i 2 arg minj log(j) log(1 hx|yji) ( ) i 2 arg minj j + c(x, yj) Lemma: With c(x, y) = log(1 hx|yi), and i := log(i), PIi(~ ) = {x 2 S2
0, c(x, yi) + i c(x, yj) + j
8j}. P1
PI3(~ )
P2 ⇢i : x 2 S2
i 1hx|yii
Wang ’04 Caffarelli-Oliker ’94
I An optimal transport problem
5
µ = probability measure on X ⌫ = prob. measure on finite Y y with density, X = manifold = P
y2Y ⌫yy
5
µ = probability measure on X ⌫ = prob. measure on finite Y T 1(y) y with density, X = manifold = P
y2Y ⌫yy
Transport map: T : X ! Y s.t. 8y 2 Y, µ(T 1({y})) = ⌫y in short: T#µ = ⌫.
5
µ = probability measure on X ⌫ = prob. measure on finite Y T 1(y) y with density, X = manifold = P
y2Y ⌫yy
Transport map: T : X ! Y s.t. 8y 2 Y, µ(T 1({y})) = ⌫y in short: T#µ = ⌫. Cc(T) = R
X c(x, T(x)) d µ(x)
Cost function: c : X ⇥ Y ! R
5
µ = probability measure on X ⌫ = prob. measure on finite Y T 1(y) y with density, X = manifold = P
y2Y ⌫yy
Transport map: T : X ! Y s.t. 8y 2 Y, µ(T 1({y})) = ⌫y in short: T#µ = ⌫. Cc(T) = R
X c(x, T(x)) d µ(x)
Cost function: c : X ⇥ Y ! R = P
y
R
T −1(y) c(x, y) d µ(x)
5
µ = probability measure on X ⌫ = prob. measure on finite Y Monge problem: Tc(µ, ⌫) := min{Cc(T); T#µ = ⌫} T 1(y) y with density, X = manifold = P
y2Y ⌫yy
Transport map: T : X ! Y s.t. 8y 2 Y, µ(T 1({y})) = ⌫y in short: T#µ = ⌫. Cc(T) = R
X c(x, T(x)) d µ(x)
Cost function: c : X ⇥ Y ! R = P
y
R
T −1(y) c(x, y) d µ(x)
Aurenhammer, Hoffman, Aronov ’98 Merigot ’2010
6
y We assume (Twist), i.e. c 2 C1 and 8x 2 X the map y 2 Y 7! rxc(x, y) is injective.
Y finite set, : Y ! R
6
y T
c (x) = arg miny2Y c(x, y) + (y)
We assume (Twist), i.e. c 2 C1 and 8x 2 X the map y 2 Y 7! rxc(x, y) is injective.
Y finite set, : Y ! R
6
y Vor
c (y) = {x 2 Rd; T c (x) = y}
T
c (x) = arg miny2Y c(x, y) + (y)
We assume (Twist), i.e. c 2 C1 and 8x 2 X the map y 2 Y 7! rxc(x, y) is injective. = generalized weighted Voronoi cell
Y finite set, : Y ! R
6
y Vor
c (y) = {x 2 Rd; T c (x) = y}
T
c (x) = arg miny2Y c(x, y) + (y)
We assume (Twist), i.e. c 2 C1 and 8x 2 X the map y 2 Y 7! rxc(x, y) is injective. NB: Under (Twist), (Vor
c (y))y2Y partitions X and T c well-defined a.e.
= generalized weighted Voronoi cell
Y finite set, : Y ! R
6
y Vor
c (y) = {x 2 Rd; T c (x) = y}
T
c (x) = arg miny2Y c(x, y) + (y)
Lemma: Given a measure µ with density and : Y ! R, the map T
c is a c-optimal transport between µ and T c#µ.
We assume (Twist), i.e. c 2 C1 and 8x 2 X the map y 2 Y 7! rxc(x, y) is injective. NB: Under (Twist), (Vor
c (y))y2Y partitions X and T c well-defined a.e.
= generalized weighted Voronoi cell
Y finite set, : Y ! R
6
y Vor
c (y) = {x 2 Rd; T c (x) = y}
T
c (x) = arg miny2Y c(x, y) + (y)
Lemma: Given a measure µ with density and : Y ! R, the map T
c is a c-optimal transport between µ and T c#µ.
We assume (Twist), i.e. c 2 C1 and 8x 2 X the map y 2 Y 7! rxc(x, y) is injective. NB: Under (Twist), (Vor
c (y))y2Y partitions X and T c well-defined a.e.
I Note: T
c#µ = P y2Y µ(Vor c (y))y.
= generalized weighted Voronoi cell
Y finite set, : Y ! R
6
y Vor
c (y) = {x 2 Rd; T c (x) = y}
T
c (x) = arg miny2Y c(x, y) + (y)
Lemma: Given a measure µ with density and : Y ! R, the map T
c is a c-optimal transport between µ and T c#µ.
We assume (Twist), i.e. c 2 C1 and 8x 2 X the map y 2 Y 7! rxc(x, y) is injective. NB: Under (Twist), (Vor
c (y))y2Y partitions X and T c well-defined a.e.
I Converse ? I Note: T
c#µ = P y2Y µ(Vor c (y))y.
= generalized weighted Voronoi cell
Y finite set, : Y ! R
7
Lemma: With c(x, y) = log(1 hx|yi), and i := log(i), PIi(~ ) = {x 2 S2
0, c(x, yi) + i c(x, yj) + j
8j}. P1
PI3(~ )
P2 y1 y2 y3
7
Lemma: With c(x, y) = log(1 hx|yi), and i := log(i), PIi(~ ) = {x 2 S2
0, c(x, yi) + i c(x, yj) + j
8j}. P1
PI3(~ )
P2 Optimal transport formulation y1 y2 y3
I PIi(~ ) = Vor
c (yi).
I T
c (x) = arg miny2Y c(x, y) + (y)
7
Lemma: With c(x, y) = log(1 hx|yi), and i := log(i), PIi(~ ) = {x 2 S2
0, c(x, yi) + i c(x, yj) + j
8j}. P1
PI3(~ )
P2 Optimal transport formulation The map T
c is a c-optimal transport between µ and T c#µ.
y1 y2 y3
I PIi(~ ) = Vor
c (yi).
I T
c (x) = arg miny2Y c(x, y) + (y)
7
Lemma: With c(x, y) = log(1 hx|yi), and i := log(i), PIi(~ ) = {x 2 S2
0, c(x, yi) + i c(x, yj) + j
8j}. P1
PI3(~ )
P2 Optimal transport formulation The map T
c is a c-optimal transport between µ and T c#µ.
y1 y2 y3
I PIi(~ ) = Vor
c (yi).
I T
c (x) = arg miny2Y c(x, y) + (y)
Problem (FF): Find 1, . . . , N such that T
c#µ = ⌫.
8
Cafarelli-Kochengin-Oliker’99: coordinate-wise ascent, with minimum increment
8
8y 2 Y \ {y0}, µ(Vorψ
c (p)) ⌫y +
Initialization: Fix y0 2 Y , let = "/N and compute s.t. While 9y 6= y0 such that µ(Vorψ
c (y)) ⌫y , do:
decrease (y) s.t. µ(Vorψ
c (y)) 2 [⌫y, ⌫y + ],
Result: s.t. for all y, |µ(Vorψ
c (y)) ⌫y| ".
Cafarelli-Kochengin-Oliker’99: coordinate-wise ascent, with minimum increment
8
8y 2 Y \ {y0}, µ(Vorψ
c (p)) ⌫y +
Initialization: Fix y0 2 Y , let = "/N and compute s.t. While 9y 6= y0 such that µ(Vorψ
c (y)) ⌫y , do:
decrease (y) s.t. µ(Vorψ
c (y)) 2 [⌫y, ⌫y + ],
Result: s.t. for all y, |µ(Vorψ
c (y)) ⌫y| ".
Cafarelli-Kochengin-Oliker’99: coordinate-wise ascent, with minimum increment
8
8y 2 Y \ {y0}, µ(Vorψ
c (p)) ⌫y +
Initialization: Fix y0 2 Y , let = "/N and compute s.t. While 9y 6= y0 such that µ(Vorψ
c (y)) ⌫y , do:
decrease (y) s.t. µ(Vorψ
c (y)) 2 [⌫y, ⌫y + ],
Result: s.t. for all y, |µ(Vorψ
c (y)) ⌫y| ".
I Complexity of SP: N 2/" steps Cafarelli-Kochengin-Oliker’99: coordinate-wise ascent, with minimum increment
8
8y 2 Y \ {y0}, µ(Vorψ
c (p)) ⌫y +
Initialization: Fix y0 2 Y , let = "/N and compute s.t. While 9y 6= y0 such that µ(Vorψ
c (y)) ⌫y , do:
decrease (y) s.t. µ(Vorψ
c (y)) 2 [⌫y, ⌫y + ],
Result: s.t. for all y, |µ(Vorψ
c (y)) ⌫y| ".
I Complexity of SP: N 2/" steps Cafarelli-Kochengin-Oliker’99: coordinate-wise ascent, with minimum increment I Generalization of Oliker–Prussner in R2 with c(x, y) = kx yk2
8
8y 2 Y \ {y0}, µ(Vorψ
c (p)) ⌫y +
Initialization: Fix y0 2 Y , let = "/N and compute s.t. While 9y 6= y0 such that µ(Vorψ
c (y)) ⌫y , do:
decrease (y) s.t. µ(Vorψ
c (y)) 2 [⌫y, ⌫y + ],
Result: s.t. for all y, |µ(Vorψ
c (y)) ⌫y| ".
I Complexity of SP: N 2/" steps I Generalization: MTW+ costs
Kitagawa ’12
Cafarelli-Kochengin-Oliker’99: coordinate-wise ascent, with minimum increment I Generalization of Oliker–Prussner in R2 with c(x, y) = kx yk2
9
Aurenhammer, Hoffman, Aronov ’98
Theorem: ~ solves (FF) iff ~ = log(~ ) maximizes Φ( ) := P
i
R
Vorψ
c (yi)[c(x, yi) + i] d µ(x) P
i i⌫i
with c(x, y) = log(1 hx|yi).
9
Aurenhammer, Hoffman, Aronov ’98
I A consequence of Kantorovich duality. Theorem: ~ solves (FF) iff ~ = log(~ ) maximizes Φ( ) := P
i
R
Vorψ
c (yi)[c(x, yi) + i] d µ(x) P
i i⌫i
with c(x, y) = log(1 hx|yi).
10
10
I @+Φ() = {v 2 Rd, Φ(µ) Φ() + hµ |vi 8µ 2 Rd}.
10
I In this case : @+Φ() = {rΦ()} a.e. I @+Φ() = {v 2 Rd, Φ(µ) Φ() + hµ |vi 8µ 2 Rd}. I Φ concave , 8 2 Rd @+Φ() 6= ;. I maximum of Φ , 0 2 @+Φ()
10
Φ( ) := P
i
R
Vorψ
c (yi)[c(x, yi) + i] d µ(x) P
i i⌫i
10
Φ( ) := P
i
R
Vorψ
c (yi)[c(x, yi) + i] d µ(x) P
i i⌫i
= R
Sd−1 min1iN[c(x, yi) + i] d µ(x) P i i⌫i
10
Φ( ) := P
i
R
Vorψ
c (yi)[c(x, yi) + i] d µ(x) P
i i⌫i
= R
Sd−1 min1iN[c(x, yi) + i] d µ(x) P i i⌫i
For all ' 2 Rd min1iN[c(x, yi) + 'i] [c(x, yTψ(x)) + 'Tψ(x)]
10
Φ( ) := P
i
R
Vorψ
c (yi)[c(x, yi) + i] d µ(x) P
i i⌫i
= R
Sd−1 min1iN[c(x, yi) + i] d µ(x) P i i⌫i
For all ' 2 Rd min1iN[c(x, yi) + 'i] [c(x, yTψ(x)) + 'Tψ(x)] T (x) = i , x 2 Vor
c (yi)
10
Φ( ) := P
i
R
Vorψ
c (yi)[c(x, yi) + i] d µ(x) P
i i⌫i
= R
Sd−1 min1iN[c(x, yi) + i] d µ(x) P i i⌫i
For all ' 2 Rd min1iN[c(x, yi) + 'i] [c(x, yTψ(x)) + 'Tψ(x)] [c(x, yTψ(x)) + Tψ(x)] + 'Tψ(x) Tψ(x) T (x) = i , x 2 Vor
c (yi)
10
Φ( ) := P
i
R
Vorψ
c (yi)[c(x, yi) + i] d µ(x) P
i i⌫i
= R
Sd−1 min1iN[c(x, yi) + i] d µ(x) P i i⌫i
For all ' 2 Rd min1iN[c(x, yi) + 'i] [c(x, yTψ(x)) + 'Tψ(x)] [c(x, yTψ(x)) + Tψ(x)] + 'Tψ(x) Tψ(x) Φ( ) + P
i i⌫i
Φ(') + P
i 'i⌫i
R
Sd−1 'Tψ(x) Tψ(x) d µ(x)
R
Sd−1
T (x) = i , x 2 Vor
c (yi)
10
Φ( ) := P
i
R
Vorψ
c (yi)[c(x, yi) + i] d µ(x) P
i i⌫i
= R
Sd−1 min1iN[c(x, yi) + i] d µ(x) P i i⌫i
Φ( ) Φ(') R
Sd−1 'Tψ(x) Tψ(x) d µ(x) P i('i i)⌫i
T (x) = i , x 2 Vor
c (yi)
10
Φ( ) := P
i
R
Vorψ
c (yi)[c(x, yi) + i] d µ(x) P
i i⌫i
= R
Sd−1 min1iN[c(x, yi) + i] d µ(x) P i i⌫i
Φ( ) Φ(') R
Sd−1 'Tψ(x) Tψ(x) d µ(x) P i('i i)⌫i
X
1iN
"Z
Vorψ
c (yi)
d µ(x) ⌫i # ('i i) T (x) = i , x 2 Vor
c (yi)
10
Φ( ) := P
i
R
Vorψ
c (yi)[c(x, yi) + i] d µ(x) P
i i⌫i
= R
Sd−1 min1iN[c(x, yi) + i] d µ(x) P i i⌫i
Φ( ) Φ(') R
Sd−1 'Tψ(x) Tψ(x) d µ(x) P i('i i)⌫i
X
1iN
"Z
Vorψ
c (yi)
d µ(x) ⌫i # ('i i) T (x) = i , x 2 Vor
c (yi)
= hDΦ( )|' i with DΦ( ) = ⇣ µ(Vor
c (yi)) ⌫i
⌘
10
Φ( ) := P
i
R
Vorψ
c (yi)[c(x, yi) + i] d µ(x) P
i i⌫i
= R
Sd−1 min1iN[c(x, yi) + i] d µ(x) P i i⌫i
Φ( ) Φ(') R
Sd−1 'Tψ(x) Tψ(x) d µ(x) P i('i i)⌫i
X
1iN
"Z
Vorψ
c (yi)
d µ(x) ⌫i # ('i i) T (x) = i , x 2 Vor
c (yi)
= hDΦ( )|' i with DΦ( ) = ⇣ µ(Vor
c (yi)) ⌫i
⌘ I DΦ( ) depends continuously on ) Φ of class C1. I DΦ( ) 2 @+Φ() ) Φ concave. I maximum of Φ , µ(Vor
c (yi)) = ⌫i 8i
11
12
Computation of descent direction / time step
LBFGS: low-storage version of the BFGS quasi-Newton scheme
I Quasi-Newton scheme:
12
R
Vorψ
c (p) d µ(x)
R
Vorψ
c (y) c(x, y) d µ(x)
Computation of descent direction / time step
LBFGS: low-storage version of the BFGS quasi-Newton scheme
Main difficulty: computation of Vor
c (y)
I Evaluation of Φ and rΦ: I Quasi-Newton scheme:
12
R
Vorψ
c (p) d µ(x)
R
Vorψ
c (y) c(x, y) d µ(x)
Computation of descent direction / time step
LBFGS: low-storage version of the BFGS quasi-Newton scheme
Main difficulty: computation of Vor
c (y)
I Evaluation of Φ and rΦ: I Quasi-Newton scheme:
13
Definition: Given P = {pi}1iN ✓ Rd and (!i)1iN 2 RN Pow!
P (pi) := {x 2 Rd; i = arg minj kx pjk2 + !j}
13
Definition: Given P = {pi}1iN ✓ Rd and (!i)1iN 2 RN Pow!
P (pi) := {x 2 Rd; i = arg minj kx pjk2 + !j}
I Efficient computation of (Pow!
P (pi))i using CGAL (d = 2, 3)
13
Lemma: With ~ = log(~ ), pi := yj
2j and !i := k yj 2j k2 1 j ,
Definition: Given P = {pi}1iN ✓ Rd and (!i)1iN 2 RN Pow!
P (pi) := {x 2 Rd; i = arg minj kx pjk2 + !j}
Vor
c (yi) = Pow! P (pi) \ S2
I Efficient computation of (Pow!
P (pi))i using CGAL (d = 2, 3)
13
Proof: x 2 Vor
c (yi) ✓ S2
= log(~ ), pi := yj
2j and !i := k yj 2j k2 1 j ,
Definition: Given P = {pi}1iN ✓ Rd and (!i)1iN 2 RN Pow!
P (pi) := {x 2 Rd; i = arg minj kx pjk2 + !j}
Vor
c (yi) = Pow! P (pi) \ S2
( ) i 2 arg minj
j 1hx|yji
I Efficient computation of (Pow!
P (pi))i using CGAL (d = 2, 3)
13
Proof: x 2 Vor
c (yi) ✓ S2
= log(~ ), pi := yj
2j and !i := k yj 2j k2 1 j ,
( ) i 2 arg minjhx| yj
j i 1 j
Definition: Given P = {pi}1iN ✓ Rd and (!i)1iN 2 RN Pow!
P (pi) := {x 2 Rd; i = arg minj kx pjk2 + !j}
Vor
c (yi) = Pow! P (pi) \ S2
( ) i 2 arg minj
j 1hx|yji
I Efficient computation of (Pow!
P (pi))i using CGAL (d = 2, 3)
13
Proof: x 2 Vor
c (yi) ✓ S2
= log(~ ), pi := yj
2j and !i := k yj 2j k2 1 j ,
( ) i 2 arg minjhx| yj
j i 1 j
( ) i 2 arg minj kx +
yj 2j k2 k yj 2j k2 1 j
pj !j
Definition: Given P = {pi}1iN ✓ Rd and (!i)1iN 2 RN Pow!
P (pi) := {x 2 Rd; i = arg minj kx pjk2 + !j}
Vor
c (yi) = Pow! P (pi) \ S2
( ) i 2 arg minj
j 1hx|yji
I Efficient computation of (Pow!
P (pi))i using CGAL (d = 2, 3)
13
Proof: x 2 Vor
c (yi) ✓ S2
= log(~ ), pi := yj
2j and !i := k yj 2j k2 1 j ,
( ) i 2 arg minjhx| yj
j i 1 j
( ) i 2 arg minj kx +
yj 2j k2 k yj 2j k2 1 j
pj !j
( ) x 2 Pow!
P (pi) \ S2
Definition: Given P = {pi}1iN ✓ Rd and (!i)1iN 2 RN Pow!
P (pi) := {x 2 Rd; i = arg minj kx pjk2 + !j}
Vor
c (yi) = Pow! P (pi) \ S2
( ) i 2 arg minj
j 1hx|yji
I Efficient computation of (Pow!
P (pi))i using CGAL (d = 2, 3)
14
I in general, the cells Ci := Pow!
P (pi) \ S2 can
be disconnected, have holes, etc.
14
I in general, the cells Ci := Pow!
P (pi) \ S2 can
be disconnected, have holes, etc. I boundary representation: a family of oriented cycles composed of circular arcs per cell.
14
I in general, the cells Ci := Pow!
P (pi) \ S2 can
be disconnected, have holes, etc. I boundary representation: a family of oriented cycles composed of circular arcs per cell. I lower complexity bound: Ω(N log N).
14
Algorithm: for each cell Ci = Pow!
P (pi) \ S2
I in general, the cells Ci := Pow!
P (pi) \ S2 can
be disconnected, have holes, etc. I boundary representation: a family of oriented cycles composed of circular arcs per cell. I lower complexity bound: Ω(N log N).
14
Algorithm: for each cell Ci = Pow!
P (pi) \ S2
every edge of Ci and S2. Set vertices V := { }. I in general, the cells Ci := Pow!
P (pi) \ S2 can
be disconnected, have holes, etc. I boundary representation: a family of oriented cycles composed of circular arcs per cell. I lower complexity bound: Ω(N log N).
14
Algorithm: for each cell Ci = Pow!
P (pi) \ S2
every edge of Ci and S2. Set vertices V := { }.
and construct oriented edges E between vertices. I in general, the cells Ci := Pow!
P (pi) \ S2 can
be disconnected, have holes, etc. I boundary representation: a family of oriented cycles composed of circular arcs per cell. I lower complexity bound: Ω(N log N).
14
Algorithm: for each cell Ci = Pow!
P (pi) \ S2
every edge of Ci and S2. Set vertices V := { }.
and construct oriented edges E between vertices.
I in general, the cells Ci := Pow!
P (pi) \ S2 can
be disconnected, have holes, etc. I boundary representation: a family of oriented cycles composed of circular arcs per cell. I lower complexity bound: Ω(N log N).
14
Algorithm: for each cell Ci = Pow!
P (pi) \ S2
every edge of Ci and S2. Set vertices V := { }.
and construct oriented edges E between vertices.
I in general, the cells Ci := Pow!
P (pi) \ S2 can
be disconnected, have holes, etc. I boundary representation: a family of oriented cycles composed of circular arcs per cell. I lower complexity bound: Ω(N log N).
14
Algorithm: for each cell Ci = Pow!
P (pi) \ S2
every edge of Ci and S2. Set vertices V := { }.
and construct oriented edges E between vertices.
I in general, the cells Ci := Pow!
P (pi) \ S2 can
be disconnected, have holes, etc. I boundary representation: a family of oriented cycles composed of circular arcs per cell. I lower complexity bound: Ω(N log N). Complexity: O(N log N + C) where C = complexity of the Power diagram.
15
⌫ = PN
i=1 ⌫ixi obtained by discretizing a picture of G. Monge.
µ = uniform measure on half-sphere S2
+
N = 1000 drawing of (Vorψ
c (yi)) (on S2 +) for = 0
15
⌫ = PN
i=1 ⌫ixi obtained by discretizing a picture of G. Monge.
µ = uniform measure on half-sphere S2
+
N = 1000 drawing of (Vorψ
c (yi)) (on S2 +) for sol
15
⌫ = PN
i=1 ⌫ixi obtained by discretizing a picture of G. Monge.
µ = uniform measure on half-sphere S2
+
N = 1000 rendering of the image reflected at infinity (using LuxRender)
16
⌫ = PN
i=1 ⌫ixi obtained by discretizing a picture of G. Monge.
µ = uniform measure on half-sphere S2
+
N = 15000 drawing of (Vorψ
c (yi)) (on S2 +) for sol
16
⌫ = PN
i=1 ⌫ixi obtained by discretizing a picture of G. Monge.
µ = uniform measure on half-sphere S2
+
N = 15000 solution to the far-field reflector problem: R(sol)
16
⌫ = PN
i=1 ⌫ixi obtained by discretizing a picture of G. Monge.
µ = uniform measure on half-sphere S2
+
N = 15000 rendering of the image reflected at infinity (using LuxRender)
17
18
Theorem: For N paraboloids, the complexity of the diagram (PIi(~ ))1iN is O(N).
18
Theorem: For N paraboloids, the complexity of the diagram (PIi(~ ))1iN is O(N). Complexity: E + F + V , where E = # edges V = # vertices F = total # of connected components
18
Theorem: For N paraboloids, the complexity of the diagram (PIi(~ ))1iN is O(N). Proof:
I F N
18
Theorem: For N paraboloids, the complexity of the diagram (PIi(~ ))1iN is O(N). Proof:
R(~ ) \ @P3(3) PI3(~ )
P1 P3
I F N
Pi Pj
Lemma: The projection of @Pi\ Pj onto the plane {y?
i } is a disc.
{y3}⊥
18
Theorem: For N paraboloids, the complexity of the diagram (PIi(~ ))1iN is O(N). Proof:
R(~ ) \ @P3(3) PI3(~ )
P1 P3
I F N
Pi Pj
Lemma: The projection of @Pi\ Pj onto the plane {y?
i } is a disc.
{y3}⊥
= ) the projection of R(~ ) \ @Pi on {yi}? is convex
18
Theorem: For N paraboloids, the complexity of the diagram (PIi(~ ))1iN is O(N). Proof:
R(~ ) \ @P3(3) PI3(~ )
P1 P3
I F N
Pi Pj
Lemma: The projection of @Pi\ Pj onto the plane {y?
i } is a disc.
{y3}⊥
= ) the projection of R(~ ) \ @Pi on {yi}? is convex = ) PIi(~ ) is connected.
18
Theorem: For N paraboloids, the complexity of the diagram (PIi(~ ))1iN is O(N). Proof:
I F N
I Every vertex has 3 edges ) 3V 2E.
18
Theorem: For N paraboloids, the complexity of the diagram (PIi(~ ))1iN is O(N). Proof:
I F N
I Every vertex has 3 edges ) 3V 2E. I Euler’s formula V E + F = 2 implies V 2F 4 and E 3F 6.
19
Proposition: Computing (PIi(~ ))i requires at least Ω(N log N) operations.
19
Proposition: Computing (PIi(~ ))i requires at least Ω(N log N) operations. Proof: reduction to a sorting problem
19
Proposition: Computing (PIi(~ ))i requires at least Ω(N log N) operations. Proof: reduction to a sorting problem
I Let t1, . . . , tN 2 R ti
19
Proposition: Computing (PIi(~ ))i requires at least Ω(N log N) operations. Proof: reduction to a sorting problem
I Let t1, . . . , tN 2 R I yi = '(ti) 2 S2 and i = cste. ti yi '
19
Proposition: Computing (PIi(~ ))i requires at least Ω(N log N) operations. Proof: reduction to a sorting problem
I Let t1, . . . , tN 2 R I yi = '(ti) 2 S2 and i = cste. I PIi(~ ) = Powω
P (yi) \ S2
with pi = yi and !i = cste. ti yi ' pi
19
Proposition: Computing (PIi(~ ))i requires at least Ω(N log N) operations. Proof: reduction to a sorting problem
I Let t1, . . . , tN 2 R I yi = '(ti) 2 S2 and i = cste. I PIi(~ ) = Powω
P (yi) \ S2
with pi = yi and !i = cste. ti yi ' pi
19
Proposition: Computing (PIi(~ ))i requires at least Ω(N log N) operations. Proof: reduction to a sorting problem
I Let t1, . . . , tN 2 R I yi = '(ti) 2 S2 and i = cste. I PIi(~ ) = Powω
P (yi) \ S2
with pi = yi and !i = cste. ti yi ' pi
I There exists a cycle of size N in the dual
I There exists a cycle of size N in the dual
19
Proposition: Computing (PIi(~ ))i requires at least Ω(N log N) operations. Proof: reduction to a sorting problem
I Let t1, . . . , tN 2 R I yi = '(ti) 2 S2 and i = cste. I PIi(~ ) = Powω
P (yi) \ S2
with pi = yi and !i = cste. ti yi ' pi
I There exists a cycle of size N in the dual
I There exists a cycle of size N in the dual
I Finding the cycle , sorting t1, . . . tN
20
21
y1 y2 y3 Punctual light at origin o, µ measure on S2
i ⌫iyi on S2 1
21
y1 y2 y3 Punctual light at origin o, µ measure on S2
i ⌫iyi on S2 1
P1
P2 µ Pi(i) = convex hull of paraboloid with focal o, direction yi and focal distance i
21
y1 y2 y3 Punctual light at origin o, µ measure on S2
i ⌫iyi on S2 1
P1
P2 µ Pi(i) = convex hull of paraboloid with focal o, direction yi and focal distance i R(~ ) = @
i=1Pi(i)
) \ @P3(3)
21
y1 y2 y3 Punctual light at origin o, µ measure on S2
i ⌫iyi on S2 1
Pi(i) = convex hull of paraboloid with focal o, direction yi and focal distance i R(~ ) = @
i=1Pi(i)
) \ @P3(3)
21
y1 y2 y3 Punctual light at origin o, µ measure on S2
i ⌫iyi on S2 1
Pi(i) = convex hull of paraboloid with focal o, direction yi and focal distance i R(~ ) = @
i=1Pi(i)
) = ⇡S2
) \ @Pi(i))
PU3(~ ) R(~ ) \ @P3(3)
21
y1 y2 y3 Punctual light at origin o, µ measure on S2
i ⌫iyi on S2 1
Pi(i) = convex hull of paraboloid with focal o, direction yi and focal distance i R(~ ) = @
i=1Pi(i)
) = ⇡S2
) \ @Pi(i)) Far-field reflector antenna problem:
PU3(~ ) R(~ ) \ @P3(3) Problem (FF’): Find 1, . . . , N such that for every i, µ(PUi(~ )) = ⌫i.
22
with pi :=
yj 2i and !i := k yi 2i k2 + 1 i .
Lemma: PUi(~ ) = Pow!
P (pi) \ S2
y1 y2 y3
PU3(~ ) R(~ ) \ @P3(3)
22
with pi :=
yj 2i and !i := k yi 2i k2 + 1 i .
Lemma: PUi(~ ) = Pow!
P (pi) \ S2
y1 y2 y3
PU3(~ ) R(~ ) \ @P3(3)
Complexity:
22
with pi :=
yj 2i and !i := k yi 2i k2 + 1 i .
Lemma: PUi(~ ) = Pow!
P (pi) \ S2
y1 y2 y3
PU3(~ ) R(~ ) \ @P3(3)
Pi Pj
I Projection of @Pi\Pj onto y?
i is a disc.
Complexity:
22
with pi :=
yj 2i and !i := k yi 2i k2 + 1 i .
Lemma: PUi(~ ) = Pow!
P (pi) \ S2
y1 y2 y3
PU3(~ ) R(~ ) \ @P3(3)
Pi Pj
I Projection of @Pi\Pj onto y?
i is a disc.
I However, the projection of R(~ ) \ @Pi is not always connected. Complexity:
22
with pi :=
yj 2i and !i := k yi 2i k2 + 1 i .
Lemma: PUi(~ ) = Pow!
P (pi) \ S2
y1 y2 y3
PU3(~ ) R(~ ) \ @P3(3)
Pi Pj
I Projection of @Pi\Pj onto y?
i is a disc.
I However, the projection of R(~ ) \ @Pi is not always connected. Complexity: ) unknown complexity
23
y1 y2 y3 Punctual light at origin o, µ measure on S2
i ⌫iyi on R3
23
y1 y2 y3 Punctual light at origin o, µ measure on S2
i ⌫iyi on R3
and yi, and eccentricity ei
23
y1 y2 y3 Punctual light at origin o, µ measure on S2
i ⌫iyi on R3
µ
and yi, and eccentricity ei
23
y1 y2 y3 Punctual light at origin o, µ measure on S2
i ⌫iyi on R3
µ
and yi, and eccentricity ei R(~ e) = @
i=1Ei(ei)
23
y1 y2 y3 Punctual light at origin o, µ measure on S2
i ⌫iyi on R3
and yi, and eccentricity ei R(~ e) = @
i=1Ei(ei)
e) = ⇡S2
e) \ @Ei(i)) EI1(~ e)
23
y1 y2 y3 Punctual light at origin o, µ measure on S2
i ⌫iyi on R3
and yi, and eccentricity ei R(~ e) = @
i=1Ei(ei)
e) = ⇡S2
e) \ @Ei(i)) Near-field reflector antenna problem: Problem (NF): Find e1, . . . , eN such that for every i, µ(EIi(~ e)) = ⌫i .
amount of light reflected to the point yi.
EI1(~ e)
Oliker ’04
24
⇢i : x 2 S2
di 1eihx|yi/kyiki where di = kyik(1e2
i )
2ei
. I @Ei(ei) is parameterized in radial coordinates by
24
⇢i : x 2 S2
di 1eihx|yi/kyiki where di = kyik(1e2
i )
2ei
. I @Ei(ei) is parameterized in radial coordinates by
e) ( ) i 2 arg minj
dj 1ejhx|yj/kyjki
24
( ) i 2 arg minj log(dj) log(1 ejhx|yj/kyjki) ⇢i : x 2 S2
di 1eihx|yi/kyiki where di = kyik(1e2
i )
2ei
. I @Ei(ei) is parameterized in radial coordinates by
e) ( ) i 2 arg minj
dj 1ejhx|yj/kyjki
24
( ) i 2 arg minj log(dj) log(1 ejhx|yj/kyjki) ⇢i : x 2 S2
di 1eihx|yi/kyiki where di = kyik(1e2
i )
2ei
. non-linearity I @Ei(ei) is parameterized in radial coordinates by
e) ( ) i 2 arg minj
dj 1ejhx|yj/kyjki
24
( ) i 2 arg minj log(dj) log(1 ejhx|yj/kyjki) ⇢i : x 2 S2
di 1eihx|yi/kyiki where di = kyik(1e2
i )
2ei
. non-linearity I @Ei(ei) is parameterized in radial coordinates by
I x 2 EIi(~ e) ( ) i 2 arg minj
dj 1ejhx|yj/kyjki
25
y1 y2
e)
25
y1 y2
e) with pi :=
eiyi 2dikyik and !i := e2
i
2d2
i 1
di .
Lemma: EIi(~ e) = Pow!
P (pi) \ S2
25
y1 y2
e) with pi :=
eiyi 2dikyik and !i := e2
i
2d2
i 1
di .
Lemma: EIi(~ e) = Pow!
P (pi) \ S2
) yi =
4 (!i+kpik2)24kpik2 pi and ei = 2kpik !ikpik2 .
I Inversion of pi :=
eiyi 2dikyik and !i :=
⇣
eiyi 2di
⌘2 1
di :
25
y1 y2
e) with pi :=
eiyi 2dikyik and !i := e2
i
2d2
i 1
di .
Lemma: EIi(~ e) = Pow!
P (pi) \ S2
) yi =
4 (!i+kpik2)24kpik2 pi and ei = 2kpik !ikpik2 .
I Inversion of pi :=
eiyi 2dikyik and !i :=
⇣
eiyi 2di
⌘2 1
di :
I ei 2 (0, 1) ( ) !i < kpik2 and !i < 1 (kpik + 1)2. Can always be satisfied by adding a negative constant to !i
25
y1 y2
e) with pi :=
eiyi 2dikyik and !i := e2
i
2d2
i 1
di .
Lemma: EIi(~ e) = Pow!
P (pi) \ S2
) yi =
4 (!i+kpik2)24kpik2 pi and ei = 2kpik !ikpik2 .
I Inversion of pi :=
eiyi 2dikyik and !i :=
⇣
eiyi 2di
⌘2 1
di :
I ei 2 (0, 1) ( ) !i < kpik2 and !i < 1 (kpik + 1)2. Can always be satisfied by adding a negative constant to !i I We can associate a family of ellipsoid to any (pi, !i)i
26
Proposition: The complexity of (Pow!
P (pi) \ S2) is Ω(N 2).
26
Proposition: The complexity of (Pow!
P (pi) \ S2) is Ω(N 2).
Proof: construct (pi)1iN s.t. the diagram Vor(pi) \ S2 has N 2 edges.
26
radius 2 " length ⌘ Proposition: The complexity of (Pow!
P (pi) \ S2) is Ω(N 2).
Proof: construct (pi)1iN s.t. the diagram Vor(pi) \ S2 has N 2 edges.
26
radius 2 " length ⌘ Proposition: The complexity of (Pow!
P (pi) \ S2) is Ω(N 2).
Proof: construct (pi)1iN s.t. the diagram Vor(pi) \ S2 has N 2 edges. I N points: k = bN/2c p1, . . . , pk : pk+1, . . . , pN :
26
radius 2 " length ⌘ Proposition: The complexity of (Pow!
P (pi) \ S2) is Ω(N 2).
Proof: construct (pi)1iN s.t. the diagram Vor(pi) \ S2 has N 2 edges. I N points: k = bN/2c I 8i, j, 9 two distinct points in Vor(pi) \ Vor(pj) \ S2 pi pj p1, . . . , pk : pk+1, . . . , pN :
26
Proposition: The complexity of (Pow!
P (pi) \ S2) is Ω(N 2).
Proof: construct (pi)1iN s.t. the diagram Vor(pi) \ S2 has N 2 edges. I N points: k = bN/2c I 8i, j, 9 two distinct points in Vor(pi) \ Vor(pj) \ S2 pi pj p1, . . . , pk : pk+1, . . . , pN :
26
Proposition: The complexity of (Pow!
P (pi) \ S2) is Ω(N 2).
Proof: construct (pi)1iN s.t. the diagram Vor(pi) \ S2 has N 2 edges. I N points: k = bN/2c I 8i, j, 9 two distinct points in Vor(pi) \ Vor(pj) \ S2 p1, . . . , pk : pk+1, . . . , pN :
26
Proposition: The complexity of (Pow!
P (pi) \ S2) is Ω(N 2).
Proof: construct (pi)1iN s.t. the diagram Vor(pi) \ S2 has N 2 edges. I N points: k = bN/2c I 8i, j, 9 two distinct points in Vor(pi) \ Vor(pj) \ S2 p1, . . . , pk : pk+1, . . . , pN : ) O(N 2) edges
26
Proposition: The complexity of (Pow!
P (pi) \ S2) is Ω(N 2).
Proof: construct (pi)1iN s.t. the diagram Vor(pi) \ S2 has N 2 edges. I N points: k = bN/2c I 8i, j, 9 two distinct points in Vor(pi) \ Vor(pj) \ S2 p1, . . . , pk : pk+1, . . . , pN : ) O(N 2) edges Similar results for the complexity of the ellipsoid union.
27
A simple quasi-Newton scheme can be used to solve rather large (15k points) geometric instances of optimal transport.
27
A simple quasi-Newton scheme can be used to solve rather large (15k points) geometric instances of optimal transport. Power diagrams can be used to compute efficiently the c-Voronoi cells
27
A simple quasi-Newton scheme can be used to solve rather large (15k points) geometric instances of optimal transport.
I quantitative stability results ? I Near field reflector problem Power diagrams can be used to compute efficiently the c-Voronoi cells I complexity of paraboloid union ? higher dimension ?
27
A simple quasi-Newton scheme can be used to solve rather large (15k points) geometric instances of optimal transport.
I quantitative stability results ? I Near field reflector problem
Power diagrams can be used to compute efficiently the c-Voronoi cells I complexity of paraboloid union ? higher dimension ?