1.5. Solving

solve is the gateway between a parsed model and everything else. You pass it a model object and a handful of options; you get back a solved model and a return code. This chapter is the reference for the option surface, the return-code discipline, and the diagnostic protocol when a solve fails.

1.5.1. The signature 

[m, retcode] = solve(m, name = value, name = value, ...);

The first input is a model object. Subsequent arguments are name-value pairs validated by an arguments block on solve.

The first output is the model object with the solution attached; downstream methods (irf, simulate, estimate, filter, …) consume it. The second output is the retcode – always capture it and always check it.

1.5.2. The retcode pattern 

The retcode is the single source of truth about whether the solve succeeded:

[m, retcode] = solve(m);
assert(retcode == 0, decipher(retcode));

retcode == 0 – the model is solved.
retcode > 0 – the solve failed in a specific, named way.
decipher(retcode) returns the human-readable explanation.

Never read the model object’s solution fields without first checking retcode. A non-zero retcode can leave the solution empty or partially populated; downstream methods will then return empty results or raise unrelated errors.

A complete list of return codes appears at the end of this chapter.

1.5.3. Picking a solver 

The solver option selects the underlying algorithm. For constant-parameter models RISE picks a sensible default. For switching models it defaults to mfi (functional iteration). You may force a specific solver by prefixing the name with +:

m = solve(m, solver = '+mn');     % force Newton, no auto-pick

The available solvers:

Name	Family	Notes
`mfi`	Functional iteration (switching default)	Cheapest per iteration. Can miss solutions where the operator is not a contraction; see the diagnostic protocol below.
`mn`	Newton (switching), no Kronecker	Robust to non-contraction. Use when `mfi` returns retcode 21.
`mnk`	Newton (switching), with Kronecker	Same robustness as `mn` with a different inner solve; sometimes faster on larger systems.
`mfi_full`, `mn_full`, `mnk_full`	Full (non-sparse) variants	Use when sparsity exploitation in the standard variants causes numerical issues.
`fwz`	Farmer-Waggoner-Zha Newton	Their algorithm; sometimes finds solutions the others miss, sometimes the reverse.
`dsge_groebner`	Groebner-basis all-solutions	Finds every solution. Can be very slow (minutes to never) on non-trivial models. Do not include in a default diagnostic sweep.
`dsge_udc`	Undetermined coefficients	For switching models.
`dsge_schur`	Schur decomposition	Only with the Maih-Waggoner perturbation strategy.
`rise_1`	Constant-parameter default	Generalised Schur / QZ decomposition.
`klein`	Klein’s method	Constant-parameter; classical reference algorithm.
`cr`, `cyclic_reduction`	Cyclic reduction	Constant-parameter.
`AIM`	AIM solver	Constant-parameter.

A user-defined solver – a function handle with a specific signature – is also accepted; see the solve doc-comment for the required inputs and outputs.

1.5.4. Solver options 

solver may be a plain string or a two-element cell {name, opts} where opts is a struct that RISE passes through to the inner iterator. The fields use MATLAB-standard option names:

opts = struct( ...
    MaxIter        = 1e6, ...
    TolFun         = 1e-10, ...
    Display        = 'iter', ...
    ObjectiveLimit = 1e+12);

[m, retcode] = solve(m, solver = {'mn', opts});

Field	Default	Role
`MaxIter`	`1000`	Maximum number of iterations.
`TolFun`	`sqrt(eps)`	Convergence tolerance on the residual.
`Display`	`'off'`	Set to `'iter'` to print convergence per iteration.
`ObjectiveLimit`	`1e+12`	Treat the iteration as divergent if the residual exceeds this.

Note

The legacy toolbox accepted the same options under the names fix_point_maxiter, fix_point_TolFun, fix_point_verbose, fix_point_explosion_limit. The modern convention is the MATLAB-standard names above; the legacy names have been retired as user-facing options and passing them now errors as an unknown option.

1.5.5. Order of approximation 

solve_order controls the order of the perturbation expansion:

m = solve(m, solve_order = 3);

Supported orders are 1 through 5. Orders above 1 require that the model was parsed with at least the matching max_deriv_order – this is a parser-time choice:

m = dsge_model('mymodel.rs', max_deriv_order = 3);
m = solve(m, solve_order = 3);

If you parsed at a lower order than you want to solve, re-parse the model – you cannot upgrade after the fact.

At order \(\geq 2\) a simulated path can diverge. Pruning is an explicit opt-in for keeping such paths bounded, off by default and enabled through the simul_pruned option; see Forecasting and simulation.

1.5.6. Optimal policy 

Two orthogonal options control optimal-policy solutions:

Option	Values	Default
`solve_policy_type`	`'ramsey'`, `'discretion'` (case-insensitive)	`'ramsey'`
`solve_policy_equilibrium`	`'OLE'`, `'MPE'` (case-insensitive)	`'OLE'`

solve_policy_type controls commitment versus discretion:

'ramsey' – full commitment; the planner re-optimizes only at t = 0.
'discretion' – time-consistent / Markov-perfect policy; the planner re-optimizes every period.

solve_policy_equilibrium controls the equilibrium concept in multi-player games:

'OLE' (open-loop equilibrium) – the classical Lagrangian formulation in which each player takes the opponent’s path as given.
'MPE' (Markov-perfect equilibrium) – each player takes the opponent’s policy function as given.

All four combinations are supported. Both options are passed to solve – they are not declared in the model file. Loose commitment and stochastic replanning are also available through a small extension of this surface; see Optimal policy in the Modeling chapter for the model-side syntax.

1.5.7. Occasionally-binding constraints 

The piecewise-linear (OccBin) solve of an occasionally-binding constraint is provided by the occbin.paradigm solve paradigm, activated through the solve_paradigm option:

C = { struct('var','R','dir',-1,'bound',1,'chain','ocb') };   % R >= 1
m = set(m,'solve_paradigm',{@occbin.paradigm,'ref',1,'constraints',C});
m = solve(m);

The model-side construction (the binding chain and the convex-combination equation) is covered in Occasionally-binding constraints, and the paradigm interface – the constraint-cell grammar and its options – in Extending RISE through paradigms.

1.5.8. Perturbation strategy 

solve_perturbation_type selects how the perturbation is taken around the deterministic steady state of a switching model. The strategies differ in which parameters and steady states are perturbed.

Value	Strategy
`'maih'` or `'m'` (default)	Maih (2014): regime-specific steady states.
`'maih_waggoner'` or `'mw'`	Maih and Waggoner (2018): regime-specific steady states plus perturbed transition probabilities on all chains.
`{'mw', mkv_list}`	As above but only the chains in `mkv_list` have their probabilities perturbed.
`'frwz'`	Foerster, Rubio-Ramirez, Waggoner and Zha (2016): unique steady state at the ergodic mean of the parameters, perturbing all switching parameters.
`{'frwz', part_list}`	As above but only the parameters in `part_list` are perturbed.

mkv_list and part_list may also be {'*none'}, {'*all'}, {'*all_but'}, or {'*all_but', 'v1', 'v2', ...}.

The default 'maih' is the right choice unless you have a specific reason to use one of the others – see the FRWZ paper for when their strategy is preferred.

Note

Under the default 'maih' perturbation, RISE auto-declares several auxiliary parameters (iota_frwz_*, ss_frwz_*, perturbator, …) that are not used in the model. Their absence from your calibration triggers a benign "1(incomplete!)" warning at parse time. This warning is not an error; do not chase it. Chase the retcode.

1.5.9. Other options worth knowing 

Option	Role
`solve_linear`	`true` to use more efficient steady-state algorithms when the model is truly linear.
`solve_bgp`	`true` to solve the model as non-stationary along a balanced growth path.
`solve_derivatives_only`	`true` to compute the symbolic derivatives but skip the solve. Useful for diagnostics.
`solve_kill_g`	`true` to kill off the growth element in the first-order solution.
`solve_automatic_differentiator`	Function handle for the automatic differentiator engine (default: `@adolm.diff`).
`solve_user_defined_shocks`	Struct describing non-Gaussian shock distributions. See the `solve` doc-comment.
`solve_qz_criterium`	Eigenvalue cutoff for constant-parameter QZ-based solvers (default: `sqrt(eps)`).

1.5.10. The diagnostic protocol when solve fails 

When solve returns a non-zero retcode, work down this list in order; each step is cheap and rules out a class of causes:

Read the retcode. decipher(retcode) gives the human-readable reason. If the message is specific (e.g. “Complex solution”, “Cannot take a log-expansion of a variable whose steady state is close to 0”), fix the model, not the solver settings.
If retcode is 21 (“Maximum Number of iterations reached or multiple solutions”), raise the iteration cap first. It is the cheapest test and rules out simple truncation:
```
opts = struct(MaxIter = 1e8);
[m, retcode] = solve(m, solver = {'mfi', opts});
```
If retcode == 0, you were just being truncated.
If retcode is still 21, switch to a Newton solver. mfi’s functional iteration is not a contraction on every model; Newton solvers are robust to that and find fixed points mfi cannot reach:
```
[m, retcode] = solve(m, solver = '+mn');
% or
[m, retcode] = solve(m, solver = '+mnk');
```
If Newton returns retcode 25 (“System unstable”), the model is genuinely unstable for this parameterization. Investigate the calibration – vary the parameters and recheck the determinacy boundary. The solver is not the issue.
If retcode is 28 (“transition matrix must be diagonal”), an OBC chain has non-zero off-diagonal transition probabilities. Set them to zero – OBC switching is endogenous, not stochastic.
If retcode is in the evaluation range (1-11), the failure is upstream of the solver. Read the message; it points at the steady-state file, the Jacobian, the transition matrix construction, or the parameter restrictions. Fix that before touching solve.

1.5.11. Return-code reference 

The canonical table lives in utils.error.registry; this is the human-readable summary, in three groups.

Evaluation (pre-solve)

Code	Meaning
`1`	Steady state could not solve
`11`	Complex or nan steady-state residuals
`2`	NaNs in Jacobian
`3`	Problem in transition matrix
`4`	Parameter restrictions violated
`5`	Definitions are NaN, Inf, or imaginary
`6`	NaNs or Infs in planner objective
`7`	Bounds or restriction violation

Solving 

Code	Meaning
`21`	Maximum number of iterations reached or multiple solutions
`210`	Discretionary outer loop: max iterations reached
`22`	NaNs in solution or no solution
`223`	Complex solution
`23`	Explosion limit reached
`230`	Discretionary outer loop: explosion limit reached
`24`	Explosive solution
`25`	System unstable
`251`	DSGE-VAR solution: system unstable
`26`	VAR approximation to the DSGE failed
`261`	DSGE prior is not proper
`262`	DSGE-VAR solution: covariance matrix not positive definite
`263`	SVAR solution: matrix not invertible
`27`	Cannot take a log-expansion of a variable whose steady state is close to 0
`28`	For this solver the transition matrix must be diagonal

Filtering and likelihood (downstream of solve, listed for completeness)

Code	Meaning
`301`	Max iterations reached in Lyapunov solution
`302`	NaNs in Lyapunov solution
`303`	Explosion limit reached in Lyapunov solution
`304`	Failure of covariance in sigma-points determination
`305`	Covariance of forecast errors not positive definite
`306`	Unlikely parameter vector
`307`	NaN/Inf/complex in log prior
`308`	Inconsistent ergodic probabilities
`309`	Error in user-defined endogenous priors
`310`	NaNs or Inf in initial conditions for the state vector
`311`	NaNs or Inf in initial covariance matrix
`312`	NaN/Inf/complex likelihood
`313`	Forecast covariance positive definite but determinant = 0

1.5.12. The retcode-aware pattern 

The pattern is the same for every solve in a script:

[m, retcode] = solve(m, solver = '+mn');

if retcode ~= 0
    warning('solve failed: %s -- trying mfi with raised MaxIter', ...
            decipher(retcode));
    opts = struct(MaxIter = 1e8);
    [m, retcode] = solve(m, solver = {'mfi', opts});
end

assert(retcode == 0, decipher(retcode));
print_solution(m);

Always capture the retcode; always decipher on failure; do not read the solution fields until retcode == 0. See Teaser Example for the canonical worked solve on the standard documentation model.