Tuning¶

Basics of tuning¶

There are only two tuning parameters: the number of internal and external iterations. We find that in many cases the tuning procedure is easy and intuitive. For problems without state constraints, only the number of internal iterations is of importance. These parameters are specified online.

At the moment, the selection of these parameters is made entirely by the user. In many embedded systems, the number of iterations may be limited by the processor computational power. More generally, the user may need to compare the MPC controller performance given by the solution of an exact solver (like CVXOPT) against that given by the solution of $\muaompc$ for a given number of iterations. For example, the comparison could be made using the stage cost at each point of a given trajectory (see [ZKF13]). In the end, the precise number of iterations strongly depends on the application.

The penalty parameter¶

An optional third tuning value is the penalty parameter $\mu$ , which is selected off-line (i.e. specified in the system module). $\muaompc$ will by default automatically compute a good value for $\mu$ if none is specified (recommended). Roughly speaking, a large penalty parameter implies that a low number of external iterations are required to reach good performance, especially when constraints are active. However, more internal iterations are necessary, because the condition number of the internal problem increases. The opposite is also true, a small $\mu$ makes the internal problem easier to solve, especially if no state constraints are active. When the state constraint are active, however, the required number of external iterations is higher.

By now it should be clear that the selection of an appropriate value of $\mu$ (not too low, not too high) is crucial.

Although in general not recommended, $\muaompc$ allows experienced users to explicitely set a value for $\mu$ in the system module. The selection of the penalty parameter $\mu$ is easily done via the function find_penalty_parameter under the tun object. For example, using the aircraft system from the tutorial:

import muaompc
mpc = muaompc.ltidt.setup_mpc_problem('sys_aircraft')
mu = mpc.tun.find_penalty_parameter()

mu now contains a value of the penalty parameter that is not too high but not too low. This is just an initial value that may help as starting point for further fine tuning. A function that may come in handy while fine tuning $\mu$ is calc_int_cn. It calculates the condition number of the internal problem for a list of penalty parameters. Continuing with the example:

mpc.tun.calc_int_cn([0.5 * mu, 2 * mu])

will return the condition number for each penalty parameter on the list. This may help to check that the selected mu does not make the internal problem ill-conditioned. Next section presents a procedure that may help reduce the condition number of the Hessian of the QP, and therefore, also of the internal problem, allowing the use of higher values of mu and making the solver faster.

Conditioning the Hessian¶

The underlying optimization algorithm used by $\muaompc$ , greatly benefits from a low condition number (say, below 100). To achieve this, we provide the function reduce_H_cn, which basically solves a nonlinear program (NLP). A description follows.

MPCTuning.reduce_H_cn(xref, uref, stw, inw, cn_ub=100, stw0=None, inw0=None)¶

Find weighting matrices that decrease the Hessian’s condition number.

Parameters:

Parameters:	xref (numpy array) – the trajectory the states should follow. `xref` has shape `(states, points)`. uref (numpy array) – the sequence of inputs for the given state trajectory `xref`. `uref` has shape `(inputs, points)`. stw (numpy array) – the relative weights for the states given as a vector of shape (states,). This vector correspond to the diagonal of the state weighting matrix `Q`, i.e. `Q = diag(stw)`. inw (numpy array) – the relative weights for the inputs given as a vector of shape (inputs,). This vector correspond to the diagonal of the input weighting matrix `R`, i.e. `R = diag(inw)`. cn_ub (scalar) – upper bound for the condition number of the Hessian. It must be greater than 1. stw0 (None or numpy array) – initial guess of the solution for the state weights. If `None`, it will be computed as an identity matrix of appropriate shape. inw0 (None or numpy array) – initial guess of the solution for the inputs weights of shape inputs If `None`, it will be computed as an identity matrix of appropriate shape.
Returns:	a dictionary with the diagonals of the newly tuned weights. The keys are as follows: “stw” is the diagonal of the weighting matrix for the states. “inw” is the diagonal of the weighting matrix for the inputs.

xref (numpy array) – the trajectory the states should follow. xref has shape (states, points).
uref (numpy array) – the sequence of inputs for the given state trajectory xref. uref has shape (inputs, points).
stw (numpy array) – the relative weights for the states given as a vector of shape (states,). This vector correspond to the diagonal of the state weighting matrix Q, i.e. Q = diag(stw).
inw (numpy array) – the relative weights for the inputs given as a vector of shape (inputs,). This vector correspond to the diagonal of the input weighting matrix R, i.e. R = diag(inw).
cn_ub (scalar) – upper bound for the condition number of the Hessian. It must be greater than 1.
stw0 (None or numpy array) – initial guess of the solution for the state weights. If None, it will be computed as an identity matrix of appropriate shape.
inw0 (None or numpy array) – initial guess of the solution for the inputs weights of shape inputs If None, it will be computed as an identity matrix of appropriate shape.

Returns:

a dictionary with the diagonals of the newly tuned weights. The keys are as follows:

“stw” is the diagonal of the weighting matrix for the states.

“inw” is the diagonal of the weighting matrix for the inputs.

Given a reference state and input trajectory (xref, uref), the difference of a simulated trajectory for the unconstrained system (with initial condition the first point of the given trajectory (xref, uref)) should be minimized with respect to the diagonal weigthing matrices of the QP. This problem is posed as a nonlinear program (NLP).

This NLP takes a reference trajectory as input (for the states and the input sequences), the reference weighting matrices that set the desired controller performance, and an upper bound on the condition number of the QP’s Hessian. The NLP returns weighting matrices that give a similar controller performance and make $\kappa_H$ (the condition number of the Hessian) lower than the specified upper bound. As this NLP is in general nonconvex, there might be several local solutions, and finding a (good) solution or not might depend on the initial guesses for the weighting matrices. For details see [WAD11].

The procedure to find suitable weighting matrices that reduce $\kappa_H$ is as follows:

Manually tune your MPC controller, as you would normally do. At the moment, only diagonal matrices are accepted. These matrices will be used as base for the NLP.

Generate a trajectory for the states and the corresponding input sequences that is representative of the MPC application.

Optionally, select an upper bound for $\kappa_H$ .

Optionally, select an initial guess for the weighting matrices.

Repeat for different parameters if the optimization was not successful, or if you are not satisfied with the results.

Take for example the system in A more complex MPC problem. We consider a change in the altitude to represent the typical behaviour of the control system. A very brief example follows:

from numpy import eye, diag
import muaompc
mpc = muaompc.ltidt.setup_mpc_problem('sys_aircraft')
# xref, uref have been already generated via simulation
# the reference weighting matrices have been manually tuned as identity
Qref = diag(eye(5))
Rref = diag(eye(1))
r = mpc.tun.reduce_H_cn(xref, uref, Qref, Rref)
Qtun = r['stw']
Rtun = r['intw']

Example¶

A complete example of how to use this function can be found in the tutorial/ltidt/tun_h_cn.