# Benchmark OpenSQP against other solvers on CUTEst problems (nx,nc<=100)


```python
'''Benchmark OpenSQP against other solvers on CUTEst problems (nx,nc<=100)'''

import numpy as np
from modopt import CUTEstProblem, optimize, OpenSQP
import pycutest
import time
import contextlib
import io
import gc # garbage collector

algs = ['SNOPT', 'IPOPT', 'PySLSQP', 'TrustConstr', 'OpenSQP']

performance = {}
history = {}
time_loop = 1
pc_prob = None
prob    = None

# Only import required problems based on the table
from modopt.benchmarking import filter_cutest_problems
prob_list = filter_cutest_problems(num_vars=[1,100], num_cons=[0,100])
print('# of problems after filtering:', len(prob_list))

# Remove some problems that cause import issues
remove_probs = ['DMN15102LS', 'DMN15103LS', 'DMN15332LS', 'DMN15333LS', 'DMN37142LS', 'DMN37143LS', 'BLEACHNG']
prob_list = [p for p in prob_list if p not in remove_probs]
print('# of problems after removing problems:', len(prob_list))

# SKIP problems for now that cause seg faults or 
# malloc errors with any of the solvers.
# Unfortunately, these problems need to be run separately 
# for solvers that have issues with them.
# This list might be different and needs to be updated 
# for different machines / OS / solver versions.
skip_probs = ['AVION2', # seg fault with IPOPT
              'BARDNE', # seg fault with OpenSQP
              'BOX3NE', # seg fault with PySLSQP
              'CHNRSBNE', # malloc error in Trust-Constr
              'MISRA1A', # malloc error in OpenSQP
              'MISRA1BLS', # seg fault
              'MISRA1C', # seg fault in OpenSQP
              'MISRA1CLS', # seg fault in OpenSQP
              'RAT42', # seg fault in OpenSQP
              'RAT43', # seg fault in PySLSQP
              ]
prob_list = [p for p in prob_list if p not in skip_probs]
print('# of problems after skipping problems:', len(prob_list))

for i, prob_name in enumerate(prob_list):
    # Import pycutest problem
    del pc_prob
    gc.collect()
    pc_prob = pycutest.import_problem(prob_name)
    print(f'[{i}.]', 'Problem name [num_vars, num_cons]:', prob_name, f'[{pc_prob.n}]', f'[{pc_prob.m}]')
    print('='*50)

    # Create modopt problem
    del prob
    gc.collect()
    prob = CUTEstProblem(cutest_problem=pc_prob)
    maxiter=250
    
    for alg in algs:
        solver = alg
        options = {}
        if alg=='IPOPT':
            options = {'hessian_approximation': 'limited-memory', 'limited_memory_max_history': 1000, 'max_iter':maxiter, 'print_level': 0, 'tol':1e-6}
        elif alg == 'TrustConstr':
            options = {'ignore_exact_hessian': True, 'maxiter': maxiter, 'gtol': 2e-5, 'xtol': 1e-100}
        elif alg == 'SNOPT':
            options = {'Verbose': False, 'Major optimality': 1.22e-4, 'Major iterations': maxiter, 'Major feasibility': 2e-6, "Verify level": -1, 
                       "Hessian frequency": 1000, "Hessian updates": 1000}
        elif alg == 'PySLSQP':
            options = {'maxiter': maxiter, 'iprint': 0, 'acc': 1e-6}
        elif alg == 'OpenSQP':
            options = {'maxiter': maxiter, 'opt_tol': 1.22e-4, 'feas_tol': 2e-6}

        print(f'\t{alg} \n\t------------------------')
        start_time = time.time()
        try:
            for i in range(time_loop-1):
                with contextlib.redirect_stdout(io.StringIO()):
                    if solver != 'OpenSQP':
                        results = optimize(prob, solver=solver, solver_options=options, recording=False, turn_off_outputs=True)
                    else:
                        results = OpenSQP(prob, maxiter=maxiter, opt_tol=options['opt_tol'], feas_tol=options['feas_tol'], recording=False, turn_off_outputs=True).solve()

            if solver != 'OpenSQP':
                results = optimize(prob, solver=solver, solver_options=options, recording=False)
            else:
                results = OpenSQP(prob, maxiter=maxiter, opt_tol=options['opt_tol'], feas_tol=options['feas_tol'], recording=False).solve()

            opt_time = (time.time() - start_time) / time_loop

            if alg == 'SNOPT':
                success = (results['info'] == 1)
            elif alg in ['IPOPT', 'IPOPT-2']:
                with open(f"{results['out_dir']}/ipopt_output.txt", 'r') as f:
                    ipopt_output = f.read()
                success = ('EXIT: Optimal Solution Found.' in ipopt_output)
            else:
                success  = results['success']
            
            nev         = prob._callback_count
            o_evals     = prob._obj_count
            g_evals     = prob._grad_count
            h_evals     = prob._hess_count
            objective   = prob._compute_objective(results['x'])

            x           = results['x']
            bd_viol_l   = np.maximum(0, prob.x_lower - x)
            bd_viol_u   = np.maximum(0, x - prob.x_upper)

            if prob.constrained:
                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)
                feasibility = np.sum(bd_viol_l) + np.sum(bd_viol_u) + np.sum(con_viol_l) + np.sum(con_viol_u)
            else:
                feasibility = np.sum(bd_viol_l) + np.sum(bd_viol_u)

            if alg == 'OpenSQP':
                nev     = results['nfev'] + results['ngev']
                nev     = nev * 2 if prob.constrained else nev
                o_evals = results['nfev']
                g_evals = results['ngev']
                h_evals = 0
            elif alg == 'SNOPT':
                nev     = nev * 4 if prob.constrained else nev * 2
                
        except Exception as e:
            print(f'Error: {e}')
            success = False
            nev = 1e6
            o_evals = 1e6
            g_evals = 1e6
            h_evals = 1e6
            objective = 1e6
            opt_time = 1e6
            feasibility = 1e6

        performance[prob.problem_name, alg] = {'time': opt_time,
                                               'success': success,
                                               'nev': nev,
                                               'objective': objective,
                                               'feasibility': feasibility}

failed_probs = {alg: [] for alg in algs}
# Print performance
print('\nPerformance')
print('='*50)
for key, value in performance.items():
    print(f"{str(key):45}:", value)
    if value['success'] == False:
        failed_probs[key[1]].append(key[0])

from modopt.benchmarking import plot_performance_profiles
plot_performance_profiles(performance, save_figname='performance.pdf')
```