This includes duplicate names for irk schemes and other things that have been moved to Enums.

Adapt the novel step eq. form also to other DCS types.

Based on #62

Re-implement functionality to calculate the type of stationarity point in RelaxationSolver.

Based on #62

NosnocIntegrator needs its own separate options.

Currently with the new solver syntax TestMpccCrossComp is a test that really tests multiple modules. We should re-architect this.

Based on discussion in #62

We want to make the front end of mpccsol more clear and as such need to consolidate the relaxation types we provide.

New README should contain a clearer explanation of the steps that NOSNOC uses to solve an OCP/simulate a system.

Currently, we support ell_1 relaxation of some terminal constraints

1. the terminal constraints of the ocp
2. some of the FESD/time-freezing terminal constraints to meet a certain final numerical/physical time.

I would still steer the relaxation of these constraints with appropriate separate options, as we have them now.

In addition to this, it would be nice to be able to relax other constraints (e.g. using vdx internally to define according variables and constraints). In particular it would be nice to relax:

1. Path constraints
2. Dynamics constraints (maybe to separate even between the differential and algebraic part of the dynamics)

Choices for relaxing:

1. ell 2 (every constraint gets a separate slack), maybe we can have upper bounds on slacks to avoid unbounded problems
2. ell1 (every constraint gets a separate slack)
3. ell_inf (every constraint group gets one slack)

Choice for penalty parameter:

1. fixed value (separate parameter for every constraint group, by groups i mean: path, dynamics, terminal, terminal_algoritmic (nice and instructive names are a plus))
2. use rho = 1/sigma - we should decide if we can mix fixed values and homotopy penalities. one one hand it gives a lot of flexibility, on the other hand too many choices might be confusing.

Currently the results struct has lots of noise in it with possibly empty entries. Remove these and document what we store in the results.

This may be irrelevant if we use vdx.

This option was experimental in the mpcc homotopy and proved not to bring much.
I suggest removing it completely as it adds no value but requires maintenance.

Based on conversation in #62.

Re-implement the fixed active set NLP extraction for finding a polished solution.

Verify correctness of reformulations via more granular unit tests.

First implementation of new design for nosnoc. For details see https://github.com/nurkanovic/nosnoc/tree/class-org-plan/organization_proposal .

clear all
clc
close all
import casadi.*
%% model parameters
e = 0.2; u_max = 9; beta = 0.0; mu=0;
%% NOSNOC settings
problem_options = NosnocProblemOptions();
solver_options = NosnocSolverOptions(); %% Optionally call this function to have an overview of all options.

problem_options.n_s = 3;
solver_options.mpcc_mode = 'Scholtes_ineq';
solver_options.homotopy_update_rule = 'superlinear';
% problem_options.step_equilibration = 'direct_homotopy';
problem_options.use_fesd = 1;
problem_options.time_freezing = 1;

problem_options.local_speed_of_time_variable = 1;
problem_options.stagewise_clock_constraint = 0;
problem_options.nonsmooth_switching_fun = 0;
problem_options.pss_lift_step_functions = 0;

problem_options.no_initial_impacts = 0;
% solver
solver_options.opts_casadi_nlp.ipopt.max_iter = 1e5;
solver_options.opts_casadi_nlp.ipopt.tol = 1e-7;
solver_options.opts_casadi_nlp.ipopt.acceptable_tol = 1e-5;
solver_options.use_previous_solution_as_initial_guess = 0;
solver_options.N_homotopy = 6;
solver_options.print_level = 3;
solver_options.comp_tol = 1e-6;
solver_options.elastic_scholtes = 1;
%% model parameters
g=9.81;
q_target = [4;0.5];

%% solver
model = NosnocModel();
model.mu_f = mu;

problem_options.T = 4;
problem_options.N_stages = 20;
problem_options.N_finite_elements = 3;
% model equations
q = SX.sym('q',2);
v = SX.sym('v',2);
u = SX.sym('u',2);

model.lbx = [0;-0.001;0;-inf]; model.ubx = [4;inf;inf;inf];
model.x0 = [0;0.5;0;0];
model.x = [q;v]; model.u = u; model.e = e ;
model.f_c = q(2)-0;
model.f_v = [u-[0;g]-betavsqrt(v(1)^2^2+v(2)^2+1e-3)];
model.J_tangent = [1; 0];
model.dims.n_dim_contact = 2;

% Objective and constraints
model.f_q = u'u; model.f_q_T = 100v'*v;
model.g_path = u'*u-u_max^2;
model.g_terminal = q-[q_target];

mpcc = NosnocMPCC(problem_options, model);
solver = NosnocSolver(mpcc, solver_options);
[results,stats] = solver.solve();

By introducing a slack variable, MPVCs are equivalent to MPCCs. This can always be done manually, but it would be maybe useful to automate this step if the user indicates that the given problem is a mpvc (or we automatically recognize it by looking at the bounds of G and H) and transform into an mpcc.
maybe it makes sense to rename everything that is mpcc into mpec.

Based on #62

Hi,
Thanks for the nice tool and quite interesting research papers.

I am not sure how to formulate one hybrid dynamics problem in order to solve it with Nosnoc. In simple terms, it is as follows:

master; % integer, influences the speed capabilities by switching the profiles 
u; % continuous control

if master==0
    f = profile_0; % define dynamics profile 0
elseif master==1
    f = profile_1; % define dynamics profile 1 
end

% goal: achieve *maximum speed*. 
% Determine best `u` (continuous control) and `master` (profile sequence e.g. profile_0, profile_1, profile_0)

I have checked examples, but I start to think you might not be dealing with these cases. I saw that you are defining the gap function. However, in this case the gap function should be kinda determined by the algorithm and it might be changing in time as well.

Would really appreciate feedback.

Best,

We have this for the plain numerical time, to make an ell_1 relaxation of the terminal time constraint.
It would be nice to have the same also for physical time constraints (maybe even also for the stagewise constraints - but to be discussed).

I want also to test with this the ease of use of vdx, once we migrated to it.

Either remove or make robust options currently listed as experimental.

Currently time-freezing only supports the Heaviside-step reformulation. Add functionality to be able to use the Stewart reformulation as well.

Currently indexes are tracked independently at different levels and must be handled somewhat manually.
Several things need to be done to make maintenance of NOSNOC easier:

• Remove FiniteElement and ControlStage classes
• Use a smarter index tracking utility to avoid pre-definition, etc. Something similar to CasADi "structures" in casadi.tools

Currently, we have for every FE an h.
Alternatively, we can define a clock state t, and have a "control variable" h_n to determine the length of a finite element.
We simply integrate the clock state as t_{n+1} = t_n + h_n; and enforce a "terminal" constraint
t_N = T_final instead of sum_n h_n = T_final
The intermediate variables t_n help eliminate sum_n h_n = T_final, which kills the block banded structure.

As we will have a new integration and OCP class after the refactoring, it would be nice to add afterwards an mpc loop that solves a sequence of OCPs.

This will also have some options, e.g. how to initialize? To use an early homotopy iterate or the solution of the previous step?
For having the new x0, use the same integrator (e..g., just shift ocp solution) or simulate with a different nosnoc integrator (e.g. higher order, more intermediate FEs).

Currently it is looped in ControlStage and it might not even work. Moreover the current for loop will be broken in newer matlab versions. R2024a already throws a warning.

It would be useful to decouple the nosnoc specific sections of the homotopy solver from the core loop and provide it as an interface to solve generic MPCCS.