Add paper.

paper/modopt_ex_16starship_benchmark.py

@@ -0,0 +1,135 @@
+'''Optimal control of Starship landing with coarsely discretized dynamics'''
+
+import numpy as np
+import matplotlib.pyplot as plt
+import time
+
+from modopt import PySLSQP, SNOPT, IPOPT, COBYQA, TrustConstr
+from modopt.postprocessing import load_variables
+from modopt.utils.profiling import profiler
+
+from examples.ex_16_4starship_landing_jax import get_problem as get_jax_prob
+
+num_steps = 20
+sol = 16.804
+obj_error = 0.168
+prob = get_jax_prob(num_steps)
+
+history      = {}
+performance  = {}
+
+def print_stats_and_save_performance(prob, solver, results, x, f,
+                                     o_time, perf_dict, hist_dict):
+    nev       = prob._callback_count
+    o_evals   = prob._obj_count
+    g_evals   = prob._grad_count
+    h_evals   = prob._hess_count
+
+    bd_viol_l  = np.maximum(0, prob.x_lower - x)
+    bd_viol_u  = np.maximum(0, x - prob.x_upper)
+
+    con         = prob._compute_constraints(x)
+    con_viol_l  = np.maximum(0, prob.c_lower - con)
+    con_viol_u  = np.maximum(0, con - prob.c_upper)
+
+    feas = np.sum(bd_viol_l) + np.sum(bd_viol_u) + np.sum(con_viol_l) + np.sum(con_viol_u)
+
+    success = np.isclose(f, sol, atol=obj_error) and feas < 1e-6
+
+    print('\tTime:', o_time, 's')
+    print('\tSuccess:', success)
+    print('\tEvaluations:', nev)
+    print('\tObj evals:', o_evals)
+    print('\tGrad evals:', g_evals)
+    print('\tHess evals:', h_evals)
+    print('\tOptimized vars:', x)
+    print('\tOptimized obj:', f)
+
+    print('\tFeasibility:', feas)
+    perf_dict[prob.problem_name, solver] = {'time': o_time,
+                                            'success': success,
+                                            'nev': nev,
+                                            'nfev': o_evals,
+                                            'ngev': g_evals,
+                                            'nhev': h_evals,
+                                            'objective': f,
+                                            'feasibility': feas}
+
+    obj_hist = load_variables(f"{results['out_dir']}/record.hdf5", 'obj')['callback_obj']
+    hist_dict[prob.problem_name, solver] = obj_hist
+
+if __name__ == '__main__':
+
+    print('\nProblem:', prob.problem_name)
+    print('='*50)
+
+    # PySLSQP
+    alg = 'PySLSQP'
+    print(f'\t{alg} \n\t------------------------')
+    optimizer = PySLSQP(prob, solver_options={'maxiter': 200, 'acc': 1e-6, 'iprint': 0}, 
+                        recording=True)
+    start_time = time.time()
+    results = optimizer.solve()
+    opt_time = time.time() - start_time
+
+    print_stats_and_save_performance(prob, alg, results, results['x'], results['objective'],
+                                     opt_time, performance, history)
+
+    # SNOPT
+    alg = 'SNOPT'
+    print(f'\t{alg} \n\t------------------------')
+    optimizer = SNOPT(prob, solver_options={'Major iterations': 200, 'Major optimality': 1e-7, 'Verbose': False},
+                        recording=True)
+
+    start_time = time.time()
+    results = optimizer.solve()
+    opt_time = time.time() - start_time
+
+    print_stats_and_save_performance(prob, alg, results, results['x'], results['objective'],
+                                     opt_time, performance, history)
+
+    # TrustConstr
+    alg = 'TrustConstr'
+    print(f'\t{alg} \n\t------------------------')
+    optimizer = TrustConstr(prob, solver_options={'maxiter': 200, 'gtol':1e-2, 'xtol':1e-6, 'ignore_exact_hessian':True}, 
+                            recording=True)
+
+    start_time = time.time()
+    results = optimizer.solve()
+    opt_time = time.time() - start_time
+
+    print_stats_and_save_performance(prob, alg, results, results['x'], results['obj'],
+                                     opt_time, performance, history)
+
+    # IPOPT
+    alg = 'IPOPT'
+    print(f'\t{alg} \n\t------------------------')
+    optimizer = IPOPT(prob, solver_options={'max_iter': 100, 'tol': 1e-3, 'print_level': 0,'accept_after_max_steps': 10},
+                      recording=True)
+
+    start_time = time.time()
+    results = optimizer.solve()
+    opt_time = time.time() - start_time
+
+    print_stats_and_save_performance(prob, alg, results, results['x'], results['f'],
+                                     opt_time, performance, history)
+
+    algs = ['PySLSQP', 'SNOPT', 'TrustConstr', 'IPOPT']
+
+    plt.figure()
+    for alg in algs:
+        y_data = history[prob.problem_name, alg]
+        plt.semilogy(y_data, label=f"{alg} ({len(y_data)})")
+    plt.xlabel('Evaluations')
+    plt.ylabel('Objective')
+    # plt.title(f'{prob.problem_name} minimization')
+    plt.legend()
+    plt.grid()
+    plt.savefig(f"{prob.problem_name}-objective-cb.png")
+    plt.close()
+
+    # Print performance
+    print('\nPerformance')
+    print('='*50)
+    for key, value in performance.items():
+        print(f"{str(key):40}:", value)

paper/paper.bib

@@ -0,0 +1,78 @@
+
+@article{kraft1988software,
+  title={A software package for sequential quadratic programming},
+  author={Kraft, Dieter},
+  journal={Forschungsbericht- Deutsche Forschungs- und Versuchsanstalt fur Luft- und Raumfahrt},
+  year={1988}
+}
+
+@article{kraft1994algorithm,
+  title={Algorithm 733: TOMP--Fortran modules for optimal control calculations},
+  author={Kraft, Dieter},
+  journal={ACM Transactions on Mathematical Software (TOMS)},
+  volume={20},
+  number={3},
+  pages={262--281},
+  year={1994},
+  publisher={ACM New York, NY, USA}
+}
+
+@article{gill2005snopt,
+  title     = {SNOPT: An SQP algorithm for large-scale constrained optimization},
+  author    = {Gill, Philip E and Murray, Walter and Saunders, Michael A},
+  journal   = {SIAM review},
+  volume    = {47},
+  number    = {1},
+  pages     = {99--131},
+  year      = {2005},
+  publisher = {SIAM},
+  doi       = {https://doi.org/10.1137/S0036144504446096}
+}
+
+@article{virtanen2020scipy,
+  title={SciPy 1.0: fundamental algorithms for scientific computing in Python},
+  author={Virtanen, Pauli and Gommers, Ralf and Oliphant, Travis E and Haberland, Matt and Reddy, Tyler and Cournapeau, David and Burovski, Evgeni and Peterson, Pearu and Weckesser, Warren and Bright, Jonathan and others},
+  journal={Nature methods},
+  volume={17},
+  number={3},
+  pages={261--272},
+  year={2020},
+  publisher={Nature Publishing Group},
+  doi       = {https://doi.org/10.1038/s41592-019-0686-2}
+}
+
+@article{wu2020pyoptsparse,
+  title={pyOptSparse: A Python framework for large-scale constrained nonlinear optimization of sparse systems},
+  author={Wu, Neil and Kenway, Gaetan and Mader, Charles A and Jasa, John and Martins, Joaquim RRA},
+  journal={Journal of Open Source Software},
+  volume={5},
+  number={54},
+  pages={2564},
+  year={2020},
+  doi={https://doi.org/10.21105/joss.02564}
+}
+
+@article{perez2012pyopt,
+  title={pyOpt: a Python-based object-oriented framework for nonlinear constrained optimization},
+  author={Perez, Ruben E and Jansen, Peter W and Martins, Joaquim RRA},
+  journal={Structural and Multidisciplinary Optimization},
+  volume={45},
+  pages={101--118},
+  year={2012},
+  publisher={Springer},
+  doi={https://doi.org/10.1007/s00158-011-0666-3}
+}
+
+@misc{NLopt,
+    author = {Steven G. Johnson},
+    title  = {The {NLopt} nonlinear-optimization package},
+    url    = {http://github.com/stevengj/nlopt},
+    year   = {2024}
+}
+
+@misc{modOpt,
+    author = {Anugrah Jo Joshy},
+    title  = {A MODular development environment and library for OPTimization algorithms},
+    url    = {https://github.com/LSDOlab/modopt},
+    year   = {2024}
+}