Hi,
I am not sure if it is user error or really a bug. But through my testing, I found that the combination of algebraic variables, nonlinear constraints, and a globalization option other than ‘FIXED_STEP’ causes the ocp_solver to crash when calling ‘solve()’.
I am using Matlab R2021a on Windows. I have installed Acados following the instruction on webpage; therefore the compiler is MinGW provided by MathWorks. I hope the following fulfills the requirement of a minimal example.
P. S. The crashes are mostly silent. However, at one point during testing with Simulink, I received an “Access violation detected” message from the MathWorks Crash Reporter. I believe this is related to the same issue. (It crashes immediately after pressing play.)
clear all;
check_acados_requirements()
%% Model
% states
x1 = casadi.MX.sym('x1');
% derivates
x1_dot = casadi.MX.sym('x1_dot');
% control
u1 = casadi.MX.sym('u1');
% algebraic variables
z1 = casadi.MX.sym('z1');
% function
eq1 = z1 - x1 + u1 - x1_dot;
eq2 = x1 + z1;
f = vertcat(eq1, eq2);
% So x_dot = -2*x1 + u1
% for constant u and x(0) = 0 -> x(t) = u/2 (1 - exp(-2t))
% Create Acados Model
model = AcadosModel();
model.x = x1;
model.xdot = x1_dot;
model.u = u1;
model.z = z1;
model.f_impl_expr = f;
%% OCP
% Create Acados OCP
ocp = AcadosOcp();
ocp.model = model;
% cost
ocp.cost.cost_type_0 = 'NONLINEAR_LS';
ocp.cost.W_0 = 1;
ocp.cost.yref_0 = 0; % for initialization only
ocp.model.cost_y_expr_0 = ocp.model.x;
ocp.cost.cost_type = 'NONLINEAR_LS';
ocp.cost.W = 1;
ocp.cost.yref = 0;
ocp.model.cost_y_expr = ocp.model.x;
ocp.cost.cost_type_e = 'NONLINEAR_LS';
ocp.cost.W_e = 1;
ocp.cost.yref_e = 0;
ocp.model.cost_y_expr_e = ocp.model.x;
% initial state constraint
ocp.constraints.x0 = 0; % required for setting initial states
% constraints on u
% ocp.constraints.constr_type_0 = 'BGH'; % is default
% ocp.constraints.constr_type = 'BGH'; % is default
% ocp.constraints.constr_type_e = 'BGH'; % is default
% ocp.constraints.idxbu = 0; % zero-based indices even for Matlab!
% ocp.constraints.lbu = -3;
% ocp.constraints.ubu = 3;
% bounds on initial
ocp.constraints.lh_0 = -3.2;
ocp.constraints.uh_0 = 3.2;
ocp.model.con_h_expr_0 = ocp.model.u;
% bounds on path
ocp.constraints.lh = -3.2;
ocp.constraints.uh = 3.2;
ocp.model.con_h_expr = ocp.model.u;
% no bounds on u at terminal node
% define solver options
ocp.solver_options.N_horizon = 10;
ocp.solver_options.tf = ocp.solver_options.N_horizon * 0.1;
ocp.solver_options.integrator_type = 'IRK';
ocp.solver_options.globalization = 'MERIT_BACKTRACKING';
% create solver
ocp_solver = AcadosOcpSolver(ocp);
%% Test
% input
ocp_solver.set('constr_x0', 0);
ocp_solver.set('cost_y_ref', 1);
ocp_solver.set('cost_y_ref_e', 1);
% solve
ocp_solver.solve();
status = ocp_solver.get('status'); % 0 - success
% output
utraj = ocp_solver.get('u');
xtraj = ocp_solver.get('x');