Hi
I use acados with the python interface and my issue refers to updating parameters after the solver is created:
I have noticed that the parameter values set via ocp.parameter_values
before the OCP solver is created effects the solution although I updated the parameters for all stages after the solver was created and before solve()
was called.
The phenomenon can be reproduced using the test_parametric_nonlinear_constaint_h.py script
Setting initial values to p_0
ocp.parameter_values = p_0
# create solver
p_set = p_0 #np.zeros(n_param)
for i in range(N):
ocp_solver.set(i, "p", p_set)
# solve()
yield the following result graph:
Setting initial values to np.zeros(n_param
and later updating with correct values p_0
ocp.parameter_values = np.zeros(n_param)
# create solver
p_set = p_0
for i in range(N):
ocp_solver.set(i, "p", p_set)
# solve()
yield slightly different results:
Any help is much appreciated, thanks!
Attached the full code to test the issue (run in folder ./acados/examples/acados_python/tests
):
import sys
sys.path.insert(0, '../pendulum_on_cart/common')
from acados_template import AcadosOcp, AcadosOcpSolver
from pendulum_model import export_pendulum_ode_model
import numpy as np
import scipy.linalg
from utils import plot_pendulum
from casadi import SX, vertcat
def main():
ocp = AcadosOcp()
# set model
model = export_pendulum_ode_model()
ocp.model = model
x = model.x
u = model.u
Tf = 1.0
nx = x.rows()
nu = u.rows()
N = 20
# set dimensions
ocp.solver_options.N_horizon = N
# set cost
Q = 2*np.diag([1e3, 1e3, 1e-2, 1e-2])
R = 2*np.diag([1e-2])
cost_W = scipy.linalg.block_diag(Q, R)
#
n_param = 42
p = SX.sym('p', n_param)
constraint_quotient = p[0]
y_param = p[1:nx+nu+1]
ocp.model.p = p
# define cost with parametric reference
ocp.cost.cost_type = 'EXTERNAL'
ocp.cost.cost_type_e = 'EXTERNAL'
residual = y_param - vertcat(x, u)
ocp.model.cost_expr_ext_cost = residual.T @ cost_W @ residual
res_e = y_param[0:nx] - x
ocp.model.cost_expr_ext_cost_e = res_e.T @ Q @ res_e
# set constraints
Fmax = 80
ocp.constraints.lh = np.array([-Fmax])
ocp.constraints.uh = np.array([+Fmax])
ocp.model.con_h_expr = model.u / constraint_quotient
ocp.constraints.lh_0 = np.array([-Fmax])
ocp.constraints.uh_0 = np.array([+Fmax])
ocp.model.con_h_expr_0 = model.u / constraint_quotient
p_0 = np.zeros(n_param)
p_0[0] = 1.0
p_0[1] = 1.0
### CHANGE THIS TO TEST THE BUG ###
# ocp.parameter_values = p_0
ocp.parameter_values = np.zeros(n_param)
ocp.constraints.x0 = np.array([0.0, np.pi, 0.0, 0.0])
ocp.solver_options.qp_solver = 'PARTIAL_CONDENSING_HPIPM' # FULL_CONDENSING_QPOASES
ocp.solver_options.hessian_approx = 'EXACT' # GAUSS_NEWTON, EXACT
ocp.solver_options.regularize_method = 'CONVEXIFY' # GAUSS_NEWTON, EXACT
ocp.solver_options.integrator_type = 'ERK'
ocp.solver_options.print_level = 0
ocp.solver_options.nlp_solver_type = 'SQP' # SQP_RTI, SQP
# set prediction horizon
ocp.solver_options.tf = Tf
# Cython
if 0:
AcadosOcpSolver.generate(ocp, json_file='acados_ocp.json')
AcadosOcpSolver.build(ocp.code_export_directory, with_cython=True)
ocp_solver = AcadosOcpSolver.create_cython_solver('acados_ocp.json')
else:
ocp_solver = AcadosOcpSolver(ocp, json_file = 'acados_ocp.json')
print("--- update params ---")
p_set = p_0 #np.zeros(n_param)
for i in range(N):
ocp_solver.set(i, "p", p_set)
solX = np.zeros((N+1, nx))
solU = np.zeros((N, nu))
print("--- solve ocp ---")
status = ocp_solver.solve()
if status != 0:
raise Exception(f'acados returned status {status}.')
# get solution
for i in range(N):
solX[i,:] = ocp_solver.get(i, "x")
solU[i,:] = ocp_solver.get(i, "u")
solX[N,:] = ocp_solver.get(N, "x")
# ocp_solver.print_statistics()
if np.any(solU > Fmax) or np.any(solU < -Fmax):
raise Exception(f"control bounds should be respected by solution, got u: {solU}")
if not np.allclose(np.max(solU), Fmax):
raise Exception(f"control should go to bounds, got u: {solU}")
if not np.allclose(np.min(solU), -Fmax):
raise Exception(f"control should go to bounds, got u: {solU}")
plot_pendulum(np.linspace(0, Tf, N+1), Fmax, solU, solX, latexify=False)
if __name__ == "__main__":
main()