1414
1515- Nonlinear optimal control of ODEs:
1616
17- $$ g(x(t_0),x(t_f)) + \int_{t_0}^{t_f} f^0(x(t), u(t))\, \mathrm{d}t \to \min $$
17+ ``` math
18+ g(x(t_0),x(t_f)) + \int_{t_0}^{t_f} f^0(x(t), u(t))\, \mathrm{d}t \to \min
19+ ```
1820
1921subject to
2022
21- $$ \dot{x}(t) = f(x(t), u(t)),\quad t \in [t_0, t_f] $$
23+ ``` math
24+ \dot{x}(t) = f(x(t), u(t)),\quad t \in [t_0, t_f]
25+ ```
2226
2327plus boundary, control and state constraints
2428
@@ -33,37 +37,347 @@ plus boundary, control and state constraints
3337
3438- Discretising an OCP into an NLP: $h_i := t_ {i+1}-t_i$,
3539
36- $$ g(X_0,X_N) + \sum_{i=0}^{N} h_i f^0(X_i,U_i) \to \min $$
40+ ``` math
41+ g(X_0,X_N) + \sum_{i=0}^{N} h_i f^0(X_i,U_i) \to \min
42+ ```
3743
3844subject to
3945
40- $$ X_{i+1} - X_i - h_i f(X_i, U_i) = 0,\quad i = 0,\dots,N-1 $$
46+ ``` math
47+ X_{i+1} - X_i - h_i f(X_i, U_i) = 0,\quad i = 0,\dots,N-1
48+ ```
4149
4250plus other constraints on $X := (X_i)_ {i=0,N}$ and $U := (U_i)_ {i=0,N}$ such as boundary and path (state and / or control) constraints :
4351
44- $$ b(t_0, X_0, t_N, X_N) = 0 $$
52+ ``` math
53+ b(t_0, X_0, t_N, X_N) = 0
54+ ```
4555
46- $$ g(t_i, X_i, U_i) = 0,\quad i = 0,\dots,N $$
56+ ``` math
57+ g(t_i, X_i, U_i) = 0,\quad i = 0,\dots,N
58+ ```
4759
4860- SIMD parallelism ($f0$, $f$, $g$) + sparsity: Kernels for GPU ([ KernelAbstraction.jl] ( https://juliagpu.github.io/KernelAbstractions.jl/stable/ ) ) and sparse linear algebra ([ CUDSS.jl] ( https://github.com/exanauts/CUDSS.jl ) )
4961- Modelling and optimising for GPU: [ ExaModels.jl] ( https://exanauts.github.io/ExaModels.jl/dev/guide ) + [ MadNLP.jl] ( https://madnlp.github.io/MadNLP.jl ) , with ** built-in AD**
5062- [ Simple example, DSL] (@ref tutorial-double-integrator-energy)
5163- Compile into an ExaModel (one pass compiler, [ syntax + semantics] ( https://github.com/control-toolbox/CTParser.jl/blob/20c6be5c953587fef10b054a95f9dc8c66b90577/src/onepass.jl#L145 ) )
52- - Simple example, generated code
64+
65+
66+ ``` @raw html
67+ <details><summary>Simple example, generated code</summary>
68+ ```
69+
5370``` julia
54- XXXX
71+ begin
72+ #= /data/caillau/CTParser.jl/src/onepass.jl:1003 =#
73+ function (; scheme = :trapezoidal , grid_size = 200 , backend = nothing , init = (0.1 , 0.1 , 0.1 ), base_type = Float64)
74+ #= /data/caillau/CTParser.jl/src/onepass.jl:1003 =#
75+ #= /data/caillau/CTParser.jl/src/onepass.jl:1004 =#
76+ LineNumberNode (0 , " box constraints: variable" )
77+ #= /data/caillau/CTParser.jl/src/onepass.jl:1005 =#
78+ begin
79+ LineNumberNode (0 , " box constraints: state" )
80+ begin
81+ var"##235" = - Inf * ones (3 )
82+ #= /data/caillau/CTParser.jl/src/onepass.jl:461 =#
83+ var"##236" = Inf * ones (3 )
84+ end
85+ end
86+ #= /data/caillau/CTParser.jl/src/onepass.jl:1006 =#
87+ begin
88+ LineNumberNode (0 , " box constraints: control" )
89+ begin
90+ var"##237" = - Inf * ones (1 )
91+ #= /data/caillau/CTParser.jl/src/onepass.jl:512 =#
92+ var"##238" = Inf * ones (1 )
93+ end
94+ end
95+ #= /data/caillau/CTParser.jl/src/onepass.jl:1007 =#
96+ var"##230" = ExaModels. ExaCore (base_type; backend = backend)
97+ #= /data/caillau/CTParser.jl/src/onepass.jl:1008 =#
98+ begin
99+ #= /data/caillau/CTParser.jl/test/test_onepass_exa.jl:23 =#
100+ var"##232" = begin
101+ #= /data/caillau/CTParser.jl/src/onepass.jl:111 =#
102+ local ex
103+ #= /data/caillau/CTParser.jl/src/onepass.jl:112 =#
104+ try
105+ #= /data/caillau/CTParser.jl/src/onepass.jl:113 =#
106+ (1 - 0 ) / grid_size
107+ catch ex
108+ #= /data/caillau/CTParser.jl/src/onepass.jl:115 =#
109+ println (" Line " , 1 , " : " , " (t ∈ [0, 1], time)" )
110+ #= /data/caillau/CTParser.jl/src/onepass.jl:116 =#
111+ throw (ex)
112+ end
113+ end
114+ #= /data/caillau/CTParser.jl/test/test_onepass_exa.jl:24 =#
115+ x = begin
116+ #= /data/caillau/CTParser.jl/src/onepass.jl:111 =#
117+ local ex
118+ #= /data/caillau/CTParser.jl/src/onepass.jl:112 =#
119+ try
120+ #= /data/caillau/CTParser.jl/src/onepass.jl:113 =#
121+ ExaModels. variable (var"##230" , 3 , 0 : grid_size; lvar = [var"##235" [i] for (i, j) = Base. product (1 : 3 , 0 : grid_size)], uvar = [var"##236" [i] for (i, j) = Base. product (1 : 3 , 0 : grid_size)], start = init[2 ])
122+ catch ex
123+ #= /data/caillau/CTParser.jl/src/onepass.jl:115 =#
124+ println (" Line " , 2 , " : " , " (x ∈ R ^ 3, state)" )
125+ #= /data/caillau/CTParser.jl/src/onepass.jl:116 =#
126+ throw (ex)
127+ end
128+ end
129+ #= /data/caillau/CTParser.jl/test/test_onepass_exa.jl:25 =#
130+ var"u##239" = begin
131+ #= /data/caillau/CTParser.jl/src/onepass.jl:111 =#
132+ local ex
133+ #= /data/caillau/CTParser.jl/src/onepass.jl:112 =#
134+ try
135+ #= /data/caillau/CTParser.jl/src/onepass.jl:113 =#
136+ ExaModels. variable (var"##230" , 1 , 0 : grid_size; lvar = [var"##237" [i] for (i, j) = Base. product (1 : 1 , 0 : grid_size)], uvar = [var"##238" [i] for (i, j) = Base. product (1 : 1 , 0 : grid_size)], start = init[3 ])
137+ catch ex
138+ #= /data/caillau/CTParser.jl/src/onepass.jl:115 =#
139+ println (" Line " , 3 , " : " , " (u ∈ R, control)" )
140+ #= /data/caillau/CTParser.jl/src/onepass.jl:116 =#
141+ throw (ex)
142+ end
143+ end
144+ #= /data/caillau/CTParser.jl/test/test_onepass_exa.jl:26 =#
145+ begin
146+ #= /data/caillau/CTParser.jl/src/onepass.jl:111 =#
147+ local ex
148+ #= /data/caillau/CTParser.jl/src/onepass.jl:112 =#
149+ try
150+ #= /data/caillau/CTParser.jl/src/onepass.jl:113 =#
151+ ExaModels. constraint (var"##230" , (x[i, 0 ] for i = 1 : 3 ); lcon = [- 1 , 0 , 0 ], ucon = [- 1 , 0 , 0 ])
152+ catch ex
153+ #= /data/caillau/CTParser.jl/src/onepass.jl:115 =#
154+ println (" Line " , 4 , " : " , " x(0) == [-1, 0, 0]" )
155+ #= /data/caillau/CTParser.jl/src/onepass.jl:116 =#
156+ throw (ex)
157+ end
158+ end
159+ #= /data/caillau/CTParser.jl/test/test_onepass_exa.jl:27 =#
160+ begin
161+ #= /data/caillau/CTParser.jl/src/onepass.jl:111 =#
162+ local ex
163+ #= /data/caillau/CTParser.jl/src/onepass.jl:112 =#
164+ try
165+ #= /data/caillau/CTParser.jl/src/onepass.jl:113 =#
166+ ExaModels. constraint (var"##230" , (x[i, grid_size] for i = 1 : 2 ); lcon = [0 , 0 ], ucon = [0 , 0 ])
167+ catch ex
168+ #= /data/caillau/CTParser.jl/src/onepass.jl:115 =#
169+ println (" Line " , 5 , " : " , " (x[1:2])(1) == [0, 0]" )
170+ #= /data/caillau/CTParser.jl/src/onepass.jl:116 =#
171+ throw (ex)
172+ end
173+ end
174+ #= /data/caillau/CTParser.jl/test/test_onepass_exa.jl:28 =#
175+ begin
176+ #= /data/caillau/CTParser.jl/src/onepass.jl:111 =#
177+ local ex
178+ #= /data/caillau/CTParser.jl/src/onepass.jl:112 =#
179+ try
180+ #= /data/caillau/CTParser.jl/src/onepass.jl:113 =#
181+ begin
182+ #= /data/caillau/CTParser.jl/src/onepass.jl:735 =#
183+ if scheme == :trapezoidal
184+ #= /data/caillau/CTParser.jl/src/onepass.jl:736 =#
185+ ExaModels. constraint (var"##230" , ((x[1 , j + 1 ] - x[1 , j]) - (var"##232" * (x[2 , j] + x[2 , j + 1 ])) / 2 for j = 0 : grid_size - 1 ))
186+ elseif #= /data/caillau/CTParser.jl/src/onepass.jl:737 =# scheme == :euler
187+ #= /data/caillau/CTParser.jl/src/onepass.jl:738 =#
188+ ExaModels. constraint (var"##230" , ((x[1 , j + 1 ] - x[1 , j]) - var"##232" * x[2 , j] for j = 0 : grid_size - 1 ))
189+ elseif #= /data/caillau/CTParser.jl/src/onepass.jl:739 =# scheme == :euler_b
190+ #= /data/caillau/CTParser.jl/src/onepass.jl:740 =#
191+ ExaModels. constraint (var"##230" , ((x[1 , j + 1 ] - x[1 , j]) - var"##232" * x[2 , j + 1 ] for j = 0 : grid_size - 1 ))
192+ else
193+ #= /data/caillau/CTParser.jl/src/onepass.jl:742 =#
194+ throw (" unknown numerical scheme" )
195+ end
196+ end
197+ catch ex
198+ #= /data/caillau/CTParser.jl/src/onepass.jl:115 =#
199+ println (" Line " , 6 , " : " , " (∂(x₁))(t) == x₂(t)" )
200+ #= /data/caillau/CTParser.jl/src/onepass.jl:116 =#
201+ throw (ex)
202+ end
203+ end
204+ #= /data/caillau/CTParser.jl/test/test_onepass_exa.jl:29 =#
205+ begin
206+ #= /data/caillau/CTParser.jl/src/onepass.jl:111 =#
207+ local ex
208+ #= /data/caillau/CTParser.jl/src/onepass.jl:112 =#
209+ try
210+ #= /data/caillau/CTParser.jl/src/onepass.jl:113 =#
211+ begin
212+ #= /data/caillau/CTParser.jl/src/onepass.jl:735 =#
213+ if scheme == :trapezoidal
214+ #= /data/caillau/CTParser.jl/src/onepass.jl:736 =#
215+ ExaModels. constraint (var"##230" , ((x[2 , j + 1 ] - x[2 , j]) - (var"##232" * (var"u##239" [1 , j] + var"u##239" [1 , j + 1 ])) / 2 for j = 0 : grid_size - 1 ))
216+ elseif #= /data/caillau/CTParser.jl/src/onepass.jl:737 =# scheme == :euler
217+ #= /data/caillau/CTParser.jl/src/onepass.jl:738 =#
218+ ExaModels. constraint (var"##230" , ((x[2 , j + 1 ] - x[2 , j]) - var"##232" * var"u##239" [1 , j] for j = 0 : grid_size - 1 ))
219+ elseif #= /data/caillau/CTParser.jl/src/onepass.jl:739 =# scheme == :euler_b
220+ #= /data/caillau/CTParser.jl/src/onepass.jl:740 =#
221+ ExaModels. constraint (var"##230" , ((x[2 , j + 1 ] - x[2 , j]) - var"##232" * var"u##239" [1 , j + 1 ] for j = 0 : grid_size - 1 ))
222+ else
223+ #= /data/caillau/CTParser.jl/src/onepass.jl:742 =#
224+ throw (" unknown numerical scheme" )
225+ end
226+ end
227+ catch ex
228+ #= /data/caillau/CTParser.jl/src/onepass.jl:115 =#
229+ println (" Line " , 7 , " : " , " (∂(x₂))(t) == u(t)" )
230+ #= /data/caillau/CTParser.jl/src/onepass.jl:116 =#
231+ throw (ex)
232+ end
233+ end
234+ #= /data/caillau/CTParser.jl/test/test_onepass_exa.jl:30 =#
235+ begin
236+ #= /data/caillau/CTParser.jl/src/onepass.jl:111 =#
237+ local ex
238+ #= /data/caillau/CTParser.jl/src/onepass.jl:112 =#
239+ try
240+ #= /data/caillau/CTParser.jl/src/onepass.jl:113 =#
241+ begin
242+ #= /data/caillau/CTParser.jl/src/onepass.jl:735 =#
243+ if scheme == :trapezoidal
244+ #= /data/caillau/CTParser.jl/src/onepass.jl:736 =#
245+ ExaModels. constraint (var"##230" , ((x[3 , j + 1 ] - x[3 , j]) - (var"##232" * (0.5 * var"u##239" [1 , j] ^ 2 + 0.5 * var"u##239" [1 , j + 1 ] ^ 2 )) / 2 for j = 0 : grid_size - 1 ))
246+ elseif #= /data/caillau/CTParser.jl/src/onepass.jl:737 =# scheme == :euler
247+ #= /data/caillau/CTParser.jl/src/onepass.jl:738 =#
248+ ExaModels. constraint (var"##230" , ((x[3 , j + 1 ] - x[3 , j]) - var"##232" * (0.5 * var"u##239" [1 , j] ^ 2 ) for j = 0 : grid_size - 1 ))
249+ elseif #= /data/caillau/CTParser.jl/src/onepass.jl:739 =# scheme == :euler_b
250+ #= /data/caillau/CTParser.jl/src/onepass.jl:740 =#
251+ ExaModels. constraint (var"##230" , ((x[3 , j + 1 ] - x[3 , j]) - var"##232" * (0.5 * var"u##239" [1 , j + 1 ] ^ 2 ) for j = 0 : grid_size - 1 ))
252+ else
253+ #= /data/caillau/CTParser.jl/src/onepass.jl:742 =#
254+ throw (" unknown numerical scheme" )
255+ end
256+ end
257+ catch ex
258+ #= /data/caillau/CTParser.jl/src/onepass.jl:115 =#
259+ println (" Line " , 8 , " : " , " (∂(x₃))(t) == 0.5 * u(t) ^ 2" )
260+ #= /data/caillau/CTParser.jl/src/onepass.jl:116 =#
261+ throw (ex)
262+ end
263+ end
264+ #= /data/caillau/CTParser.jl/test/test_onepass_exa.jl:31 =#
265+ begin
266+ #= /data/caillau/CTParser.jl/src/onepass.jl:111 =#
267+ local ex
268+ #= /data/caillau/CTParser.jl/src/onepass.jl:112 =#
269+ try
270+ #= /data/caillau/CTParser.jl/src/onepass.jl:113 =#
271+ ExaModels. objective (var"##230" , x[3 , grid_size])
272+ catch ex
273+ #= /data/caillau/CTParser.jl/src/onepass.jl:115 =#
274+ println (" Line " , 9 , " : " , " x₃(1) → min" )
275+ #= /data/caillau/CTParser.jl/src/onepass.jl:116 =#
276+ throw (ex)
277+ end
278+ end
279+ end
280+ #= /data/caillau/CTParser.jl/src/onepass.jl:1009 =#
281+ begin
282+ #= /data/caillau/CTParser.jl/src/onepass.jl:994 =#
283+ ! (isempty ([1 , 2 , 3 ])) || throw (CTBase. ParsingError (" dynamics not defined" ))
284+ #= /data/caillau/CTParser.jl/src/onepass.jl:995 =#
285+ sort ([1 , 2 , 3 ]) == 1 : 3 || throw (CTBase. ParsingError (" some coordinates of dynamics undefined" ))
286+ end
287+ #= /data/caillau/CTParser.jl/src/onepass.jl:1010 =#
288+ return ExaModels. ExaModel (var"##230" )
289+ end
290+ end
55291```
292+
293+ ``` @raw html
294+ </details>
295+ ```
296+
56297- ** Remark.** Automated scalarisation of (linear) range constraints
57298``` julia
58299XXXX
59300```
60301- Solving (MadNLP + CUDSS)
61302``` julia
62- XXXX
303+ This is MadNLP version v0.8.7 , running with cuDSS v0.4.0
304+
305+ Number of nonzeros in constraint Jacobian............ : 12005
306+ Number of nonzeros in Lagrangian Hessian............. : 9000
307+
308+ Total number of variables............................ : 4004
309+ variables with only lower bounds: 0
310+ variables with lower and upper bounds: 0
311+ variables with only upper bounds: 0
312+ Total number of equality constraints................. : 3005
313+ Total number of inequality constraints............... : 0
314+ inequality constraints with only lower bounds: 0
315+ inequality constraints with lower and upper bounds: 0
316+ inequality constraints with only upper bounds: 0
317+
318+ iter objective inf_pr inf_du lg (mu) || d|| lg (rg) alpha_du alpha_pr ls
319+ 0 1.0000000e-01 1.10e+00 1.00e+00 - 1.0 0.00e+00 - 0.00e+00 0.00e+00 0
320+ 1 1.0001760e-01 1.10e+00 3.84e-03 - 1.0 6.88e+02 - 4.0 1.00e+00 2.00e-07 h 2
321+ 2 - 5.2365072e-03 1.89e-02 1.79e-07 - 1.0 6.16e+00 - 4.5 1.00e+00 1.00e+00 h 1
322+ 3 5.9939627e+00 2.28e-03 1.66e-04 - 3.8 6.00e+00 - 5.0 9.99e-01 1.00e+00 h 1
323+ 4 5.9996194e+00 3.00e-06 1.05e-06 - 3.8 7.77e-02 - 1.00e+00 1.00e+00 h 1
324+
325+ Number of Iterations.... : 4
326+
327+ (scaled) (unscaled)
328+ Objective............... : 5.9996194116826898e+00 5.9996194116826898e+00
329+ Dual infeasibility...... : 1.0504883830561350e-06 1.0504883830561350e-06
330+ Constraint violation.... : 2.9972943093481582e-06 2.9972943093481582e-06
331+ Complementarity......... : 2.0007458990898193e-06 2.0007458990898193e-06
332+ Overall NLP error....... : 2.9972943093481582e-06 2.9972943093481582e-06
333+
334+ Number of objective function evaluations = 6
335+ Number of objective gradient evaluations = 5
336+ Number of constraint evaluations = 6
337+ Number of constraint Jacobian evaluations = 5
338+ Number of Lagrangian Hessian evaluations = 4
339+ Total wall- clock secs in solver (w/ o fun. eval./ lin. alg.) = 0.423
340+ Total wall- clock secs in linear solver = 0.055
341+ Total wall- clock secs in NLP function evaluations = 0.005
342+ Total wall- clock secs = 0.483
63343```
64344- [ Goddard problem] ( https://control-toolbox.org/Tutorials.jl/stable/tutorial-goddard/ )
65345``` julia
66- XXXX
346+ This is MadNLP version v0.8.4 , running with cuDSS v0.3.0
347+
348+ Number of nonzeros in constraint Jacobian............ : 135017
349+ Number of nonzeros in Lagrangian Hessian............. : 130008
350+
351+ Total number of variables............................ : 35008
352+ variables with only lower bounds: 0
353+ variables with lower and upper bounds: 0
354+ variables with only upper bounds: 0
355+ Total number of equality constraints................. : 30007
356+ Total number of inequality constraints............... : 15004
357+ inequality constraints with only lower bounds: 5002
358+ inequality constraints with lower and upper bounds: 10002
359+ inequality constraints with only upper bounds: 0
360+
361+ [... ]
362+
363+ Number of Iterations.... : 35
364+
365+ (scaled) (unscaled)
366+ Objective............... : - 1.0142336978192805e+00 - 1.0142336978192805e+00
367+ Dual infeasibility...... : 4.7384318691001681e-13 4.7384318691001681e-13
368+ Constraint violation.... : 1.4068322357215250e-09 1.4068322357215250e-09
369+ Complementarity......... : 9.0909295306344959e-09 9.0909295306344959e-09
370+ Overall NLP error....... : 9.0909295306344959e-09 9.0909295306344959e-09
371+
372+ Number of objective function evaluations = 36
373+ Number of objective gradient evaluations = 36
374+ Number of constraint evaluations = 36
375+ Number of constraint Jacobian evaluations = 36
376+ Number of Lagrangian Hessian evaluations = 35
377+ Total wall- clock secs in solver (w/ o fun. eval./ lin. alg.) = 0.911
378+ Total wall- clock secs in linear solver = 0.227
379+ Total wall- clock secs in NLP function evaluations = 0.059
380+ Total wall- clock secs = 1.198
67381```
68382
69383## Wrap up
0 commit comments