ipopt_sys/
lib.rs

1//   Copyright 2018 Egor Larionov
2//
3//   Licensed under the Apache License, Version 2.0 (the "License");
4//   you may not use this file except in compliance with the License.
5//   You may obtain a copy of the License at
6//
7//       http://www.apache.org/licenses/LICENSE-2.0
8//
9//   Unless required by applicable law or agreed to in writing, software
10//   distributed under the License is distributed on an "AS IS" BASIS,
11//   WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
12//   See the License for the specific language governing permissions and
13//   limitations under the License.
14
15#![allow(non_upper_case_globals)]
16#![allow(non_camel_case_types)]
17#![allow(non_snake_case)]
18include!(concat!(env!("OUT_DIR"), "/ipopt_cnlp.rs"));
19
20#[cfg(test)]
21mod tests {
22    use super::*;
23    use approx::*;
24    use std::ffi::CString;
25    use std::slice;
26
27    /// A small structure used to store state between calls to Ipopt NLP callbacks.
28    #[derive(Debug)]
29    struct UserData {
30        g_offset: [CNLP_Number; 2],
31    }
32
33    /// Test Ipopt raw bindings. This will also serve as an example of the raw C API.
34    #[test]
35    fn hs071_test() {
36        // rough comparator
37        let approx_eq =
38            |a: f64, b: f64| assert_relative_eq!(a, b, max_relative = 1e-6, epsilon = 1e-14);
39
40        let mut nlp: CNLP_ProblemPtr = ::std::ptr::null_mut();
41        let create_status = unsafe {
42            /* create the CNLP_Problem */
43            cnlp_create_problem(
44                &mut nlp as *mut CNLP_ProblemPtr,
45                0,
46                Some(sizes),
47                Some(init),
48                Some(bounds),
49                Some(eval_f),
50                Some(eval_g),
51                Some(eval_grad_f),
52                Some(eval_jac_g),
53                Some(eval_h),
54                None,
55            )
56        };
57
58        assert_eq!(create_status, CNLP_CreateProblemStatus_CNLP_SUCCESS);
59
60        /* Set some options */
61        let tol_str = CString::new("tol").unwrap();
62        let mu_strategy_str = CString::new("mu_strategy").unwrap();
63        let adaptive_str = CString::new("adaptive").unwrap();
64        let print_lvl_str = CString::new("print_level").unwrap();
65        let sb_str = CString::new("sb").unwrap();
66        let yes_str = CString::new("yes").unwrap();
67
68        unsafe {
69            cnlp_add_int_option(nlp, print_lvl_str.as_ptr(), 0);
70            cnlp_add_num_option(nlp, tol_str.as_ptr(), 1e-7);
71            cnlp_add_str_option(nlp, mu_strategy_str.as_ptr(), adaptive_str.as_ptr());
72            cnlp_add_str_option(nlp, sb_str.as_ptr(), yes_str.as_ptr());
73            cnlp_set_intermediate_callback(nlp, Some(intermediate_cb));
74        }
75
76        let mut user_data = UserData {
77            g_offset: [0.0, 0.0],
78        };
79        let udata_ptr = (&mut user_data) as *mut UserData;
80
81        /* solve the problem */
82        let sol = unsafe {
83            cnlp_solve(
84                nlp,
85                udata_ptr as CNLP_UserDataPtr, // Pointer to user data. This will be passed unmodified
86                                               // to the callback functions.
87            ) // Problem that is to be optimized.
88        };
89
90        assert_eq!(
91            sol.status,
92            CNLP_ApplicationReturnStatus_CNLP_USER_REQUESTED_STOP
93        );
94
95        let m = 2;
96        let n = 4;
97
98        let mut g = vec![0.0; m];
99
100        // Get solution data
101        let mut x = vec![0.0; n];
102        let mut mult_g = vec![0.0; m];
103        let mut mult_x_L = vec![0.0; n];
104        let mut mult_x_U = vec![0.0; n];
105
106        x.copy_from_slice(unsafe { slice::from_raw_parts(sol.data.x, n) });
107        g.copy_from_slice(unsafe { slice::from_raw_parts(sol.g, m) });
108        mult_g.copy_from_slice(unsafe { slice::from_raw_parts(sol.data.mult_g, m) });
109        mult_x_L.copy_from_slice(unsafe { slice::from_raw_parts(sol.data.mult_x_L, n) });
110        mult_x_U.copy_from_slice(unsafe { slice::from_raw_parts(sol.data.mult_x_U, n) });
111
112        approx_eq(x[0], 1.000000e+00);
113        approx_eq(x[1], 4.743000e+00);
114        approx_eq(x[2], 3.821150e+00);
115        approx_eq(x[3], 1.379408e+00);
116
117        approx_eq(mult_g[0], -5.522936e-01);
118        approx_eq(mult_g[1], 1.614685e-01);
119
120        approx_eq(mult_x_L[0], 1.087872e+00);
121        approx_eq(mult_x_L[1], 4.635819e-09);
122        approx_eq(mult_x_L[2], 9.087447e-09);
123        approx_eq(mult_x_L[3], 8.555955e-09);
124        approx_eq(mult_x_U[0], 4.470027e-09);
125        approx_eq(mult_x_U[1], 4.075231e-07);
126        approx_eq(mult_x_U[2], 1.189791e-08);
127        approx_eq(mult_x_U[3], 6.398749e-09);
128
129        approx_eq(sol.obj_val, 1.701402e+01);
130
131        // Now we are going to solve this problem again, but with slightly modified
132        // constraints.  We change the constraint offset of the first constraint a bit,
133        // and resolve the problem using the warm start option.
134
135        user_data.g_offset[0] = 0.2;
136
137        let mut warm_start_str = CString::new("warm_start_init_point").unwrap();
138        let mut yes_str = CString::new("yes").unwrap();
139        let mut bound_push_str = CString::new("bound_push").unwrap();
140        let mut bound_frac_str = CString::new("bound_frac").unwrap();
141
142        unsafe {
143            cnlp_add_str_option(
144                nlp,
145                (&mut warm_start_str).as_ptr() as *mut i8,
146                (&mut yes_str).as_ptr() as *mut i8,
147            );
148            cnlp_add_num_option(nlp, (&mut bound_push_str).as_ptr() as *mut i8, 1e-5);
149            cnlp_add_num_option(nlp, (&mut bound_frac_str).as_ptr() as *mut i8, 1e-5);
150            cnlp_set_intermediate_callback(nlp, None);
151        }
152
153        let sol = unsafe { cnlp_solve(nlp, udata_ptr as CNLP_UserDataPtr) };
154
155        // Write solutions back to our managed Vecs
156        x.copy_from_slice(unsafe { slice::from_raw_parts(sol.data.x, n) });
157        g.copy_from_slice(unsafe { slice::from_raw_parts(sol.g, m) });
158        mult_g.copy_from_slice(unsafe { slice::from_raw_parts(sol.data.mult_g, m) });
159        mult_x_L.copy_from_slice(unsafe { slice::from_raw_parts(sol.data.mult_x_L, n) });
160        mult_x_U.copy_from_slice(unsafe { slice::from_raw_parts(sol.data.mult_x_U, n) });
161
162        assert_eq!(
163            sol.status,
164            CNLP_ApplicationReturnStatus_CNLP_SOLVE_SUCCEEDED
165        );
166
167        approx_eq(x[0], 1.000000e+00);
168        approx_eq(x[1], 4.749269e+00);
169        approx_eq(x[2], 3.817510e+00);
170        approx_eq(x[3], 1.367870e+00);
171
172        approx_eq(mult_g[0], -5.517016e-01);
173        approx_eq(mult_g[1], 1.592915e-01);
174
175        approx_eq(mult_x_L[0], 1.090362e+00);
176        approx_eq(mult_x_L[1], 2.664877e-12);
177        approx_eq(mult_x_L[2], 3.556758e-12);
178        approx_eq(mult_x_L[3], 2.718148e-11);
179        approx_eq(mult_x_U[0], 2.498100e-12);
180        approx_eq(mult_x_U[1], 4.013263e-11);
181        approx_eq(mult_x_U[2], 8.455630e-12);
182        approx_eq(mult_x_U[3], 2.755724e-12);
183
184        approx_eq(sol.obj_val, 1.690362e+01);
185
186        /* free allocated memory */
187        unsafe {
188            cnlp_free_problem(nlp);
189        }
190    }
191
192    /* Function Implementations */
193    unsafe extern "C" fn sizes(
194        n: *mut CNLP_Index,
195        m: *mut CNLP_Index,
196        nnz_jac_g: *mut CNLP_Index,
197        nnz_h_lag: *mut CNLP_Index,
198        _user_data: CNLP_UserDataPtr,
199    ) -> CNLP_Bool {
200        *n = 4; // CNLP_Number of variables
201        *m = 2; // CNLP_Number of constraints
202        *nnz_jac_g = 8; // CNLP_Number of jacobian non-zeros
203        *nnz_h_lag = 10; // CNLP_Number of hessian non-zeros
204        true as CNLP_Bool
205    }
206
207    unsafe extern "C" fn init(
208        n: CNLP_Index,
209        init_x: CNLP_Bool,
210        x: *mut CNLP_Number,
211        init_z: CNLP_Bool,
212        z_L: *mut CNLP_Number,
213        z_U: *mut CNLP_Number,
214        m: CNLP_Index,
215        init_lambda: CNLP_Bool,
216        lambda: *mut CNLP_Number,
217        _user_data: CNLP_UserDataPtr,
218    ) -> CNLP_Bool {
219        assert_eq!(n, 4);
220        assert_eq!(m, 2);
221        if init_x != 0 {
222            /* initialize values for the initial point */
223            *x.offset(0) = 1.0;
224            *x.offset(1) = 5.0;
225            *x.offset(2) = 5.0;
226            *x.offset(3) = 1.0;
227        }
228        if init_z != 0 {
229            /* initialize multipliers for the variable bounds */
230            for i in 0..4 {
231                *z_L.offset(i) = 0.0;
232                *z_U.offset(i) = 0.0;
233            }
234        }
235        if init_lambda != 0 {
236            *lambda.offset(0) = 0.0;
237            *lambda.offset(1) = 0.0;
238        }
239        true as CNLP_Bool
240    }
241
242    unsafe extern "C" fn bounds(
243        n: CNLP_Index,
244        x_l: *mut CNLP_Number,
245        x_u: *mut CNLP_Number,
246        m: CNLP_Index,
247        g_l: *mut CNLP_Number,
248        g_u: *mut CNLP_Number,
249        _user_data: CNLP_UserDataPtr,
250    ) -> CNLP_Bool {
251        assert!(n == 4);
252        assert!(m == 2);
253
254        /* Set the values of the constraint bounds */
255        *g_l.offset(0) = 25.0;
256        *g_l.offset(1) = 40.0;
257        *g_u.offset(0) = 2.0e19;
258        *g_u.offset(1) = 40.0;
259
260        /* Set the values for the variable bounds */
261        for i in 0..n as isize {
262            *x_l.offset(i) = 1.0;
263            *x_u.offset(i) = 5.0;
264        }
265
266        true as CNLP_Bool
267    }
268
269    unsafe extern "C" fn eval_f(
270        n: CNLP_Index,
271        x: *const CNLP_Number,
272        _new_x: CNLP_Bool,
273        obj_value: *mut CNLP_Number,
274        _user_data: CNLP_UserDataPtr,
275    ) -> CNLP_Bool {
276        assert!(n == 4);
277
278        *obj_value = *x.offset(0) * *x.offset(3) * (*x.offset(0) + *x.offset(1) + *x.offset(2))
279            + *x.offset(2);
280
281        true as CNLP_Bool
282    }
283
284    unsafe extern "C" fn eval_grad_f(
285        n: CNLP_Index,
286        x: *const CNLP_Number,
287        _new_x: CNLP_Bool,
288        grad_f: *mut CNLP_Number,
289        _user_data: CNLP_UserDataPtr,
290    ) -> CNLP_Bool {
291        assert!(n == 4);
292
293        *grad_f.offset(0) = *x.offset(0) * *x.offset(3)
294            + *x.offset(3) * (*x.offset(0) + *x.offset(1) + *x.offset(2));
295        *grad_f.offset(1) = *x.offset(0) * *x.offset(3);
296        *grad_f.offset(2) = *x.offset(0) * *x.offset(3) + 1.0;
297        *grad_f.offset(3) = *x.offset(0) * (*x.offset(0) + *x.offset(1) + *x.offset(2));
298
299        true as CNLP_Bool
300    }
301
302    unsafe extern "C" fn eval_g(
303        n: CNLP_Index,
304        x: *const CNLP_Number,
305        _new_x: CNLP_Bool,
306        m: CNLP_Index,
307        g: *mut CNLP_Number,
308        user_data_ptr: CNLP_UserDataPtr,
309    ) -> CNLP_Bool {
310        //struct MyUserData* my_data = user_data;
311
312        assert!(n == 4);
313        assert!(m == 2);
314
315        let user_data = &*(user_data_ptr as *mut UserData);
316        *g.offset(0) =
317            *x.offset(0) * *x.offset(1) * *x.offset(2) * *x.offset(3) + user_data.g_offset[0];
318        *g.offset(1) = *x.offset(0) * *x.offset(0)
319            + *x.offset(1) * *x.offset(1)
320            + *x.offset(2) * *x.offset(2)
321            + *x.offset(3) * *x.offset(3)
322            + user_data.g_offset[1];
323
324        true as CNLP_Bool
325    }
326
327    unsafe extern "C" fn eval_jac_g(
328        _n: CNLP_Index,
329        x: *const CNLP_Number,
330        _new_x: CNLP_Bool,
331        _m: CNLP_Index,
332        _nele_jac: CNLP_Index,
333        iRow: *mut CNLP_Index,
334        jCol: *mut CNLP_Index,
335        values: *mut CNLP_Number,
336        _user_data: CNLP_UserDataPtr,
337    ) -> CNLP_Bool {
338        if values.is_null() {
339            /* return the structure of the jacobian */
340
341            /* this particular jacobian is dense */
342            *iRow.offset(0) = 0;
343            *jCol.offset(0) = 0;
344            *iRow.offset(1) = 0;
345            *jCol.offset(1) = 1;
346            *iRow.offset(2) = 0;
347            *jCol.offset(2) = 2;
348            *iRow.offset(3) = 0;
349            *jCol.offset(3) = 3;
350            *iRow.offset(4) = 1;
351            *jCol.offset(4) = 0;
352            *iRow.offset(5) = 1;
353            *jCol.offset(5) = 1;
354            *iRow.offset(6) = 1;
355            *jCol.offset(6) = 2;
356            *iRow.offset(7) = 1;
357            *jCol.offset(7) = 3;
358        } else {
359            /* return the values of the jacobian of the constraints */
360
361            *values.offset(0) = *x.offset(1) * *x.offset(2) * *x.offset(3); /* 0,0 */
362            *values.offset(1) = *x.offset(0) * *x.offset(2) * *x.offset(3); /* 0,1 */
363            *values.offset(2) = *x.offset(0) * *x.offset(1) * *x.offset(3); /* 0,2 */
364            *values.offset(3) = *x.offset(0) * *x.offset(1) * *x.offset(2); /* 0,3 */
365
366            *values.offset(4) = 2.0 * *x.offset(0); /* 1,0 */
367            *values.offset(5) = 2.0 * *x.offset(1); /* 1,1 */
368            *values.offset(6) = 2.0 * *x.offset(2); /* 1,2 */
369            *values.offset(7) = 2.0 * *x.offset(3); /* 1,3 */
370        }
371
372        true as CNLP_Bool
373    }
374
375    unsafe extern "C" fn eval_h(
376        _n: CNLP_Index,
377        x: *const CNLP_Number,
378        _new_x: CNLP_Bool,
379        obj_factor: CNLP_Number,
380        _m: CNLP_Index,
381        lambda: *const CNLP_Number,
382        _new_lambda: CNLP_Bool,
383        nele_hess: CNLP_Index,
384        iRow: *mut CNLP_Index,
385        jCol: *mut CNLP_Index,
386        values: *mut CNLP_Number,
387        _user_data: CNLP_UserDataPtr,
388    ) -> CNLP_Bool {
389        if values.is_null() {
390            /* return the structure. This is a symmetric matrix, fill the lower left
391             * triangle only. */
392
393            /* the hessian for this problem is actually dense */
394            let mut idx = 0; /* nonzero element counter */
395            for row in 0..4 {
396                for col in 0..row + 1 {
397                    *iRow.offset(idx) = row;
398                    *jCol.offset(idx) = col;
399                    idx += 1;
400                }
401            }
402
403            assert!(idx == nele_hess as isize);
404        } else {
405            /* return the values. This is a symmetric matrix, fill the lower left
406             * triangle only */
407
408            /* fill the objective portion */
409            *values.offset(0) = obj_factor * (2.0 * *x.offset(3)); /* 0,0 */
410
411            *values.offset(1) = obj_factor * (*x.offset(3)); /* 1,0 */
412            *values.offset(2) = 0.0; /* 1,1 */
413
414            *values.offset(3) = obj_factor * (*x.offset(3)); /* 2,0 */
415            *values.offset(4) = 0.0; /* 2,1 */
416            *values.offset(5) = 0.0; /* 2,2 */
417
418            *values.offset(6) = obj_factor * (2.0 * *x.offset(0) + *x.offset(1) + *x.offset(2)); /* 3,0 */
419            *values.offset(7) = obj_factor * (*x.offset(0)); /* 3,1 */
420            *values.offset(8) = obj_factor * (*x.offset(0)); /* 3,2 */
421            *values.offset(9) = 0.0; /* 3,3 */
422
423            /* add the portion for the first constraint */
424            *values.offset(1) += *lambda.offset(0) * (*x.offset(2) * *x.offset(3)); /* 1,0 */
425
426            *values.offset(3) += *lambda.offset(0) * (*x.offset(1) * *x.offset(3)); /* 2,0 */
427            *values.offset(4) += *lambda.offset(0) * (*x.offset(0) * *x.offset(3)); /* 2,1 */
428
429            *values.offset(6) += *lambda.offset(0) * (*x.offset(1) * *x.offset(2)); /* 3,0 */
430            *values.offset(7) += *lambda.offset(0) * (*x.offset(0) * *x.offset(2)); /* 3,1 */
431            *values.offset(8) += *lambda.offset(0) * (*x.offset(0) * *x.offset(1)); /* 3,2 */
432
433            /* add the portion for the second constraint */
434            *values.offset(0) += *lambda.offset(1) * 2.0; /* 0,0 */
435
436            *values.offset(2) += *lambda.offset(1) * 2.0; /* 1,1 */
437
438            *values.offset(5) += *lambda.offset(1) * 2.0; /* 2,2 */
439
440            *values.offset(9) += *lambda.offset(1) * 2.0; /* 3,3 */
441        }
442
443        true as CNLP_Bool
444    }
445
446    extern "C" fn intermediate_cb(
447        _alg_mod: CNLP_AlgorithmMode,
448        _iter_count: CNLP_Index,
449        _obj_value: CNLP_Number,
450        inf_pr: CNLP_Number,
451        _inf_du: CNLP_Number,
452        _mu: CNLP_Number,
453        _d_norm: CNLP_Number,
454        _regularization_size: CNLP_Number,
455        _alpha_du: CNLP_Number,
456        _alpha_pr: CNLP_Number,
457        _ls_trials: CNLP_Index,
458        _user_data: CNLP_UserDataPtr,
459    ) -> CNLP_Bool {
460        (inf_pr >= 1e-4) as CNLP_Bool
461    }
462}
ipopt_sys/lib.rs

ipopt_sys/
lib.rs
