Home > matpower7.1 > lib > mpoption.m

mpoption

PURPOSE ^

MPOPTION Used to set and retrieve a MATPOWER options struct.

SYNOPSIS ^

function opt = mpoption(varargin)

DESCRIPTION ^

MPOPTION  Used to set and retrieve a MATPOWER options struct.

   OPT = MPOPTION
       Returns the default options struct.

   OPT = MPOPTION(OVERRIDES)
       Returns the default options struct, with some fields overridden
       by values from OVERRIDES, which can be a struct or the name of
       a function that returns a struct.

   OPT = MPOPTION(NAME1, VALUE1, NAME2, VALUE2, ...)
       Same as previous, except override options are specified by NAME,
       VALUE pairs. This can be used to set any part of the options
       struct. The names can be individual fields or multi-level field
       names with embedded periods. The values can be scalars or structs.

       For backward compatibility, the NAMES and VALUES may correspond
       to old-style MATPOWER option names (elements in the old-style
       options vector) as well.

   OPT = MPOPTION(OPT0)
       Converts an old-style options vector OPT0 into the corresponding
       options struct. If OPT0 is an options struct it does nothing.

   OPT = MPOPTION(OPT0, OVERRIDES)
       Applies overrides to an existing set of options, OPT0, which
       can be an old-style options vector or an options struct.

   OPT = MPOPTION(OPT0, NAME1, VALUE1, NAME2, VALUE2, ...)
       Same as above except it uses the old-style options vector OPT0
       as a base instead of the old default options vector.

   OPT_VECTOR = MPOPTION(OPT, [])
       Creates and returns an old-style options vector from an
       options struct OPT.

   Note: The use of old-style MATPOWER options vectors and their
         names and values has been deprecated and will be removed
         in a future version of MATPOWER. Until then, all uppercase
         option names are not permitted for new top-level options.

   Examples:
       mpopt = mpoption('pf.alg', 'FDXB', 'pf.tol', 1e-4);
       mpopt = mpoption(mpopt, 'opf.dc.solver', 'CPLEX', 'verbose', 2);

The currently defined options are as follows:

   name                    default     description [options]
----------------------    ---------   ----------------------------------
Model options:
   model                   'AC'        AC vs. DC power flow model
       [ 'AC' - use nonlinear AC model & corresponding algorithms/options  ]
       [ 'DC' - use linear DC model & corresponding algorithms/options     ]

Power Flow options:
   pf.alg                  'NR'        AC power flow algorithm
       [ 'NR'    - Newton's method (formulation depends on values of       ]
       [           pf.current_balance and pf.v_cartesian options)          ]
       [ 'NR-SP' - Newton's method (power mismatch, polar)                 ]
       [ 'NR-SC' - Newton's method (power mismatch, cartesian)             ]
       [ 'NR-SH' - Newton's method (power mismatch, hybrid)                ]
       [ 'NR-IP' - Newton's method (current mismatch, polar)               ]
       [ 'NR-IC' - Newton's method (current mismatch, cartesian)           ]
       [ 'NR-IH' - Newton's method (current mismatch, hybrid)              ]
       [ 'FDXB'  - Fast-Decoupled (XB version)                             ]
       [ 'FDBX'  - Fast-Decoupled (BX version)                             ]
       [ 'GS'    - Gauss-Seidel                                            ]
       [ 'PQSUM' - Power Summation method (radial networks only)           ]
       [ 'ISUM'  - Current Summation method (radial networks only)         ]
       [ 'YSUM'  - Admittance Summation method (radial networks only)      ]
   pf.current_balance     0           type of nodal balance equation
       [ 0 - use complex power balance equations                           ]
       [ 1 - use complex current balance equations                         ]
   pf.v_cartesian         0           voltage representation
       [ 0 - bus voltage variables represented in polar coordinates        ]
       [ 1 - bus voltage variables represented in cartesian coordinates    ]
       [ 2 - hybrid, polar updates computed via modified cartesian Jacobian]
   pf.tol                  1e-8        termination tolerance on per unit
                                       P & Q mismatch
   pf.nr.max_it            10          maximum number of iterations for
                                       Newton's method
   pf.nr.lin_solver        ''          linear solver passed to MPLINSOLVE to
                                       solve Newton update step
       [ ''      - default to '\' for small systems, 'LU3' for larger ones ]
       [ '\'     - built-in backslash operator                             ]
       [ 'LU'    - explicit default LU decomposition and back substitution ]
       [ 'LU3'   - 3 output arg form of LU, Gilbert-Peierls algorithm with ]
       [           approximate minimum degree (AMD) reordering             ]
       [ 'LU4'   - 4 output arg form of LU, UMFPACK solver (same as 'LU')  ]
       [ 'LU5'   - 5 output arg form of LU, UMFPACK solver, w/row scaling  ]
       [       (see MPLINSOLVE for complete list of all options)           ]
   pf.fd.max_it            30          maximum number of iterations for
                                       fast decoupled method
   pf.gs.max_it            1000        maximum number of iterations for
                                       Gauss-Seidel method
   pf.radial.max_it        20          maximum number of iterations for
                                       radial power flow methods
   pf.radial.vcorr         0           perform voltage correction procedure
                                       in distribution power flow
       [  0 - do NOT perform voltage correction                            ]
       [  1 - perform voltage correction                                   ]
   pf.enforce_q_lims       0           enforce gen reactive power limits at
                                       expense of |V|
       [  0 - do NOT enforce limits                                        ]
       [  1 - enforce limits, simultaneous bus type conversion             ]
       [  2 - enforce limits, one-at-a-time bus type conversion            ]

Continuation Power Flow options:
   cpf.parameterization    3           choice of parameterization
       [  1 - natural                                                      ]
       [  2 - arc length                                                   ]
       [  3 - pseudo arc length                                            ]
   cpf.stop_at             'NOSE'      determins stopping criterion
       [ 'NOSE'     - stop when nose point is reached                      ]
       [ 'FULL'     - trace full nose curve                                ]
       [ <lam_stop> - stop upon reaching specified target lambda value     ]
   cpf.enforce_p_lims      0           enforce gen active power limits
       [  0 - do NOT enforce limits                                        ]
       [  1 - enforce limits, simultaneous bus type conversion             ]
   cpf.enforce_q_lims      0           enforce gen reactive power limits at
                                       expense of |V|
       [  0 - do NOT enforce limits                                        ]
       [  1 - enforce limits, simultaneous bus type conversion             ]
   cpf.enforce_v_lims      0           enforce bus voltage magnitude limits
       [  0 - do NOT enforce limits                                        ]
       [  1 - enforce limits, termination on detection                     ]
   cpf.enforce_flow_lims   0           enforce branch flow MVA limits
       [  0 - do NOT enforce limits                                        ]
       [  1 - enforce limits, termination on detection                     ]
   cpf.step                0.05        continuation power flow step size
   cpf.adapt_step          0           toggle adaptive step size feature
       [  0 - adaptive step size disabled                                  ]
       [  1 - adaptive step size enabled                                   ]
   cpf.step_min            1e-4        minimum allowed step size
   cpf.step_max            0.2         maximum allowed step size
   cpf.adapt_step_damping  0.7         damping factor for adaptive step
                                       sizing
   cpf.adapt_step_tol      1e-3        tolerance for adaptive step sizing
   cpf.target_lam_tol      1e-5        tolerance for target lambda detection
   cpf.nose_tol            1e-5        tolerance for nose point detection (pu)
   cpf.p_lims_tol          0.01        tolerance for generator active
                                       power limit enforcement (MW)
   cpf.q_lims_tol          0.01        tolerance for generator reactive
                                       power limit enforcement (MVAR)
   cpf.v_lims_tol          1e-4        tolerance for bus voltage
                                       magnitude enforcement (p.u)
   cpf.flow_lims_tol       0.01        tolerance for line MVA flow
                                       enforcement (MVA)
   cpf.plot.level          0           control plotting of noze curve
       [  0 - do not plot nose curve                                       ]
       [  1 - plot when completed                                          ]
       [  2 - plot incrementally at each iteration                         ]
       [  3 - same as 2, with 'pause' at each iteration                    ]
   cpf.plot.bus            <empty>     index of bus whose voltage is to be
                                       plotted
   cpf.user_callback       <empty>     string containing the name of a user
                                       callback function, or struct with
                                       function name, and optional priority
                                       and/or args, or cell array of such
                                       strings and/or structs, see
                                       'help cpf_default_callback' for details

Optimal Power Flow options:
   name                    default     description [options]
----------------------    ---------   ----------------------------------
   opf.ac.solver           'DEFAULT'   AC optimal power flow solver
       [ 'DEFAULT' - choose default solver, i.e. 'MIPS'                    ]
       [ 'MIPS'    - MIPS, MATPOWER Interior Point Solver, primal/dual     ]
       [             interior point method (pure MATLAB/Octave)            ]
       [ 'FMINCON' - MATLAB Optimization Toolbox, FMINCON                  ]
       [ 'IPOPT'   - IPOPT, requires MEX interface to IPOPT solver         ]
       [             available from:                                       ]
       [                   https://github.com/coin-or/Ipopt                ]
       [ 'KNITRO'  - Artelys Knitro, requires Artelys Knitro solver,       ]
       [             available from:https://www.artelys.com/solvers/knitro/]
       [ 'MINOPF'  - MINOPF, MINOS-based solver, requires optional         ]
       [             MEX-based MINOPF package, available from:             ]
       [                   http://www.pserc.cornell.edu/minopf/            ]
       [ 'PDIPM'   - PDIPM, primal/dual interior point method, requires    ]
       [             optional MEX-based TSPOPF package, available from:    ]
       [                   http://www.pserc.cornell.edu/tspopf/            ]
       [ 'SDPOPF'  - SDPOPF, solver based on semidefinite relaxation of    ]
       [             OPF problem, requires optional packages:              ]
       [               SDP_PF, available in extras/sdp_pf                  ]
       [               YALMIP, available from:                             ]
       [                   https://yalmip.github.io                        ]
       [               SDP solver such as SeDuMi, available from:          ]
       [                   http://sedumi.ie.lehigh.edu/                    ]
       [ 'TRALM'   - TRALM, trust region based augmented Langrangian       ]
       [             method, requires TSPOPF (see 'PDIPM')                 ]
   opf.dc.solver           'DEFAULT'   DC optimal power flow solver
       [ 'DEFAULT' - choose solver based on availability in the following  ]
       [             order: 'GUROBI', 'CPLEX', 'MOSEK', 'OT',              ]
       [             'GLPK' (linear costs only), 'BPMPD', 'MIPS'           ]
       [ 'MIPS'    - MIPS, MATPOWER Interior Point Solver, primal/dual     ]
       [             interior point method (pure MATLAB/Octave)            ]
       [ 'BPMPD'   - BPMPD, requires optional MEX-based BPMPD_MEX package  ]
       [             available from: http://www.pserc.cornell.edu/bpmpd/   ]
       [ 'CLP'     - CLP, requires interface to COIN-OP LP solver          ]
       [             available from:https://github.com/coin-or/Clp         ]
       [ 'CPLEX'   - CPLEX, requires CPLEX solver available from:          ]
       [             https://www.ibm.com/analytics/cplex-optimizer         ]
       [ 'GLPK'    - GLPK, requires interface to GLPK solver               ]
       [             available from: https://www.gnu.org/software/glpk/    ]
       [             (GLPK does not work with quadratic cost functions)    ]
       [ 'GUROBI'  - GUROBI, requires Gurobi optimizer (v. 5+)             ]
       [             available from: https://www.gurobi.com/               ]
       [ 'IPOPT'   - IPOPT, requires MEX interface to IPOPT solver         ]
       [             available from:                                       ]
       [                 https://github.com/coin-or/Ipopt                  ]
       [ 'MOSEK'   - MOSEK, requires MATLAB interface to MOSEK solver      ]
       [             available from: https://www.mosek.com/                ]
       [ 'OSQP'    - OSQP, requires MATLAB interface to OSQP solver        ]
       [             available from: https://osqp.org/                     ]
       [ 'OT'      - MATLAB Optimization Toolbox, QUADPROG, LINPROG        ]
   opf.current_balance     0           type of nodal balance equation
       [ 0 - use complex power balance equations                           ]
       [ 1 - use complex current balance equations                         ]
   opf.v_cartesian         0           voltage representation
       [ 0 - bus voltage variables represented in polar coordinates        ]
       [ 1 - bus voltage variables represented in cartesian coordinates    ]
   opf.violation           5e-6        constraint violation tolerance
   opf.use_vg              0           respect gen voltage setpt     [ 0-1 ]
       [ 0 - use specified bus Vmin & Vmax, and ignore gen Vg              ]
       [ 1 - replace specified bus Vmin & Vmax by corresponding gen Vg     ]
       [ between 0 and 1 - use a weighted average of the 2 options         ]
   opf.flow_lim            'S'         quantity limited by branch flow
                                       constraints
       [ 'S' - apparent power flow (limit in MVA)                          ]
       [ 'P' - active power flow, implemented using P   (limit in MW)      ]
       [ '2' - active power flow, implemented using P^2 (limit in MW)      ]
       [ 'I' - current magnitude (limit in MVA at 1 p.u. voltage)          ]
   opf.ignore_angle_lim    0           angle diff limits for branches
       [ 0 - include angle difference limits, if specified                 ]
       [ 1 - ignore angle difference limits even if specified              ]
   opf.softlims.default    1           behavior of OPF soft limits for
                                       which parameters are not explicitly
                                       provided
       [ 0 - do not include softlims if not explicitly specified           ]
       [ 1 - include softlims w/default values if not explicitly specified ]
   opf.init_from_mpc       -1          (DEPRECATED: use opf.start instead)
                                       specify whether to use current state
                                       in MATPOWER case to initialize OPF
                                       (currently only supported by fmincon,
                                        Ipopt, Knitro and MIPS solvers,
                                        others always use mpc)
       [  -1 - MATPOWER decides, based on solver/algorithm                 ]
       [   0 - ignore current state in MATPOWER case (only applies to      ]
       [       fmincon, Ipopt, Knitro and MIPS, which use an interior pt   ]
       [       estimate; others use current state as with opf.start = 2)   ]
       [   1 - use current state in MATPOWER case                          ]
   opf.start               0           strategy for initializing OPF start pt
       [   0 - default, MATPOWER decides based on solver                   ]
       [       (currently identical to 1)                                  ]
       [   1 - ignore current state in MATPOWER case (only applies to      ]
       [       fmincon, Ipopt, Knitro and MIPS, which use an interior pt   ]
       [       estimate; others use current state as with opf.start = 2)   ]
       [   2 - use current state in MATPOWER case                          ]
       [   3 - solve power flow and use resulting state                    ]
   opf.return_raw_der      0           for AC OPF, return constraint and
                                       derivative info in results.raw
                                       (in fields g, dg, df, d2f) [ 0 or 1 ]

Output options:
   name                    default     description [options]
----------------------    ---------   ----------------------------------
   verbose                 1           amount of progress info printed
       [   0 - print no progress info                                      ]
       [   1 - print a little progress info                                ]
       [   2 - print a lot of progress info                                ]
       [   3 - print all progress info                                     ]
   out.all                 -1          controls pretty-printing of results
       [  -1 - individual flags control what prints                        ]
       [   0 - do not print anything (overrides individual flags, ignored  ]
       [       for files specified as FNAME arg to runpf(), runopf(), etc.)]
       [   1 - print everything (overrides individual flags)               ]
   out.sys_sum             1           print system summary       [ 0 or 1 ]
   out.area_sum            0           print area summaries       [ 0 or 1 ]
   out.bus                 1           print bus detail           [ 0 or 1 ]
   out.branch              1           print branch detail        [ 0 or 1 ]
   out.gen                 0           print generator detail     [ 0 or 1 ]
   out.lim.all             -1          controls constraint info output
       [  -1 - individual flags control what constraint info prints        ]
       [   0 - no constraint info (overrides individual flags)             ]
       [   1 - binding constraint info (overrides individual flags)        ]
       [   2 - all constraint info (overrides individual flags)            ]
   out.lim.v               1           control voltage limit info
       [   0 - do not print                                                ]
       [   1 - print binding constraints only                              ]
       [   2 - print all constraints                                       ]
       [   (same options for OUT_LINE_LIM, OUT_PG_LIM, OUT_QG_LIM)         ]
   out.lim.line            1           control line flow limit info
   out.lim.pg              1           control gen active power limit info
   out.lim.qg              1           control gen reactive pwr limit info
   out.force               0           print results even if success
                                       flag = 0                   [ 0 or 1 ]
   out.suppress_detail     -1          suppress all output but system summary
       [  -1 - suppress details for large systems (> 500 buses)            ]
       [   0 - do not suppress any output specified by other flags         ]
       [   1 - suppress all output except system summary section           ]
       [       (overrides individual flags, but not out.all = 1)           ]

Solver specific options:
       name                    default     description [options]
   -----------------------    ---------   ----------------------------------
   MIPS:
       mips.linsolver          ''          linear system solver
           [   '' or '\'   build-in backslash \ operator (e.g. x = A \ b)  ]
           [   'PARDISO'   PARDISO solver (if available)                   ]
       mips.feastol            0           feasibility (equality) tolerance
                                           (set to opf.violation by default)
       mips.gradtol            1e-6        gradient tolerance
       mips.comptol            1e-6        complementary condition
                                           (inequality) tolerance
       mips.costtol            1e-6        optimality tolerance
       mips.max_it             150         maximum number of iterations
       mips.step_control       0           enable step-size cntrl [ 0 or 1 ]
       mips.sc.red_it          20          maximum number of reductions per
                                           iteration with step control
       mips.xi                 0.99995     constant used in alpha updates*
       mips.sigma              0.1         centering parameter*
       mips.z0                 1           used to initialize slack variables*
       mips.alpha_min          1e-8        returns "Numerically Failed" if
                                           either alpha parameter becomes
                                           smaller than this value*
       mips.rho_min            0.95        lower bound on rho_t*
       mips.rho_max            1.05        upper bound on rho_t*
       mips.mu_threshold       1e-5        KT multipliers smaller than this
                                           value for non-binding constraints
                                           are forced to zero
       mips.max_stepsize       1e10        returns "Numerically Failed" if the
                                           2-norm of the reduced Newton step
                                           exceeds this value*
           * See the corresponding Appendix in the manual for details.

   CPLEX:
       cplex.lpmethod          0           solution algorithm for LP problems
           [   0 - automatic: let CPLEX choose                             ]
           [   1 - primal simplex                                          ]
           [   2 - dual simplex                                            ]
           [   3 - network simplex                                         ]
           [   4 - barrier                                                 ]
           [   5 - sifting                                                 ]
           [   6 - concurrent (dual, barrier, and primal)                  ]
       cplex.qpmethod          0           solution algorithm for QP problems
           [   0 - automatic: let CPLEX choose                             ]
           [   1 - primal simplex optimizer                                ]
           [   2 - dual simplex optimizer                                  ]
           [   3 - network optimizer                                       ]
           [   4 - barrier optimizer                                       ]
       cplex.opts              <empty>     see CPLEX_OPTIONS for details
       cplex.opt_fname         <empty>     see CPLEX_OPTIONS for details
       cplex.opt               0           see CPLEX_OPTIONS for details

   FMINCON:
       fmincon.alg             4           algorithm used by fmincon() for OPF
                                           for Opt Toolbox 4 and later
            [  1 - active-set (not suitable for large problems)            ]
            [  2 - interior-point, w/default 'bfgs' Hessian approx         ]
            [  3 - interior-point, w/ 'lbfgs' Hessian approx               ]
            [  4 - interior-point, w/exact user-supplied Hessian           ]
            [  5 - interior-point, w/Hessian via finite differences        ]
            [  6 - sqp (not suitable for large problems)                   ]
       fmincon.tol_x           1e-4        termination tol on x
       fmincon.tol_f           1e-4        termination tol on f
       fmincon.max_it          0           maximum number of iterations
                                                           [  0 => default ]

   GUROBI:
       gurobi.method           0           solution algorithm (Method)
           [  -1 - automatic, let Gurobi decide                            ]
           [   0 - primal simplex                                          ]
           [   1 - dual simplex                                            ]
           [   2 - barrier                                                 ]
           [   3 - concurrent (LP only)                                    ]
           [   4 - deterministic concurrent (LP only)                      ]
       gurobi.timelimit        Inf         maximum time allowed (TimeLimit)
       gurobi.threads          0           max number of threads (Threads)
       gurobi.opts             <empty>     see GUROBI_OPTIONS for details
       gurobi.opt_fname        <empty>     see GUROBI_OPTIONS for details
       gurobi.opt              0           see GUROBI_OPTIONS for details

   IPOPT:
       ipopt.opts              <empty>     see IPOPT_OPTIONS for details
       ipopt.opt_fname         <empty>     see IPOPT_OPTIONS for details
       ipopt.opt               0           see IPOPT_OPTIONS for details

   KNITRO:
       knitro.tol_x            1e-4        termination tol on x
       knitro.tol_f            1e-4        termination tol on f
       knitro.maxit            0           maximum number of iterations
                                                           [  0 => default ]
       knitro.opt_fname        <empty>     name of user-supplied native
                                           Knitro options file that overrides
                                           all other options
       knitro.opt              0           if knitro.opt_fname is empty and
                                           knitro.opt is a non-zero integer N
                                           then knitro.opt_fname is auto-
                                           generated as:
                                           'knitro_user_options_N.txt'

   LINPROG:
       linprog                 <empty>     LINPROG options passed to
                                           OPTIMOPTIONS or OPTIMSET.
                                           see LINPROG in the Optimization
                                           Toolbox for details

   MINOPF:
       minopf.feastol          0 (1e-3)    primal feasibility tolerance
                                           (set to opf.violation by default)
       minopf.rowtol           0 (1e-3)    row tolerance
       minopf.xtol             0 (1e-4)    x tolerance
       minopf.majdamp          0 (0.5)     major damping parameter
       minopf.mindamp          0 (2.0)     minor damping parameter
       minopf.penalty          0 (1.0)     penalty parameter
       minopf.major_it         0 (200)     major iterations
       minopf.minor_it         0 (2500)    minor iterations
       minopf.max_it           0 (2500)    iterations limit
       minopf.verbosity        -1          amount of progress info printed
           [  -1 - controlled by 'verbose' option                          ]
           [   0 - print nothing                                           ]
           [   1 - print only termination status message                   ]
           [   2 - print termination status and screen progress            ]
           [   3 - print screen progress, report file (usually fort.9)     ]
       minopf.core             0 (1200*nb + 2*(nb+ng)^2) memory allocation
       minopf.supbasic_lim     0 (2*nb + 2*ng) superbasics limit
       minopf.mult_price       0 (30)      multiple price

   MOSEK:
       mosek.lp_alg            0           solution algorithm
                                               (MSK_IPAR_OPTIMIZER)
           for MOSEK 8.x ...        (see MOSEK_SYMBCON for a "better way")
           [   0 - automatic: let MOSEK choose                             ]
           [   1 - dual simplex                                            ]
           [   2 - automatic: let MOSEK choose                             ]
           [   3 - automatic simplex (MOSEK chooses which simplex method)  ]
           [   4 - interior point                                          ]
           [   6 - primal simplex                                          ]
       mosek.max_it            0 (400)     interior point max iterations
                                               (MSK_IPAR_INTPNT_MAX_ITERATIONS)
       mosek.gap_tol           0 (1e-8)    interior point relative gap tol
                                               (MSK_DPAR_INTPNT_TOL_REL_GAP)
       mosek.max_time          0 (-1)      maximum time allowed
                                               (MSK_DPAR_OPTIMIZER_MAX_TIME)
       mosek.num_threads       0 (1)       max number of threads
                                               (MSK_IPAR_INTPNT_NUM_THREADS)
       mosek.opts              <empty>     see MOSEK_OPTIONS for details
       mosek.opt_fname         <empty>     see MOSEK_OPTIONS for details
       mosek.opt               0           see MOSEK_OPTIONS for details

   OSQP:
       osqp.opts               <empty>     see OSQP_OPTIONS for details

   QUADPROG:
       quadprog                <empty>     QUADPROG options passed to
                                           OPTIMOPTIONS or OPTIMSET.
                                           see QUADPROG in the Optimization
                                           Toolbox for details

   TSPOPF:
       pdipm.feastol           0           feasibility (equality) tolerance
                                           (set to opf.violation by default)
       pdipm.gradtol           1e-6        gradient tolerance
       pdipm.comptol           1e-6        complementary condition
                                           (inequality) tolerance
       pdipm.costtol           1e-6        optimality tolerance
       pdipm.max_it            150         maximum number of iterations
       pdipm.step_control      0           enable step-size cntrl [ 0 or 1 ]
       pdipm.sc.red_it         20          maximum number of reductions per
                                           iteration with step control
       pdipm.sc.smooth_ratio   0.04        piecewise linear curve smoothing
                                           ratio

       tralm.feastol           0           feasibility tolerance
                                           (set to opf.violation by default)
       tralm.primaltol         5e-4        primal variable tolerance
       tralm.dualtol           5e-4        dual variable tolerance
       tralm.costtol           1e-5        optimality tolerance
       tralm.major_it          40          maximum number of major iterations
       tralm.minor_it          40          maximum number of minor iterations
       tralm.smooth_ratio      0.04        piecewise linear curve smoothing
                                           ratio

Experimental Options:
   exp.sys_wide_zip_loads.pw   <empty>     1 x 3 vector of active load fraction
                                           to be modeled as constant power,
                                           constant current and constant
                                           impedance, respectively, where
                                           <empty> means use [1 0 0]
   exp.sys_wide_zip_loads.qw   <empty>     same for reactive power, where
                                           <empty> means use same value as
                                           for 'pw'

CROSS-REFERENCE INFORMATION ^

This function calls: This function is called by:

SUBFUNCTIONS ^

SOURCE CODE ^

0001 function opt = mpoption(varargin)
0002 %MPOPTION  Used to set and retrieve a MATPOWER options struct.
0003 %
0004 %   OPT = MPOPTION
0005 %       Returns the default options struct.
0006 %
0007 %   OPT = MPOPTION(OVERRIDES)
0008 %       Returns the default options struct, with some fields overridden
0009 %       by values from OVERRIDES, which can be a struct or the name of
0010 %       a function that returns a struct.
0011 %
0012 %   OPT = MPOPTION(NAME1, VALUE1, NAME2, VALUE2, ...)
0013 %       Same as previous, except override options are specified by NAME,
0014 %       VALUE pairs. This can be used to set any part of the options
0015 %       struct. The names can be individual fields or multi-level field
0016 %       names with embedded periods. The values can be scalars or structs.
0017 %
0018 %       For backward compatibility, the NAMES and VALUES may correspond
0019 %       to old-style MATPOWER option names (elements in the old-style
0020 %       options vector) as well.
0021 %
0022 %   OPT = MPOPTION(OPT0)
0023 %       Converts an old-style options vector OPT0 into the corresponding
0024 %       options struct. If OPT0 is an options struct it does nothing.
0025 %
0026 %   OPT = MPOPTION(OPT0, OVERRIDES)
0027 %       Applies overrides to an existing set of options, OPT0, which
0028 %       can be an old-style options vector or an options struct.
0029 %
0030 %   OPT = MPOPTION(OPT0, NAME1, VALUE1, NAME2, VALUE2, ...)
0031 %       Same as above except it uses the old-style options vector OPT0
0032 %       as a base instead of the old default options vector.
0033 %
0034 %   OPT_VECTOR = MPOPTION(OPT, [])
0035 %       Creates and returns an old-style options vector from an
0036 %       options struct OPT.
0037 %
0038 %   Note: The use of old-style MATPOWER options vectors and their
0039 %         names and values has been deprecated and will be removed
0040 %         in a future version of MATPOWER. Until then, all uppercase
0041 %         option names are not permitted for new top-level options.
0042 %
0043 %   Examples:
0044 %       mpopt = mpoption('pf.alg', 'FDXB', 'pf.tol', 1e-4);
0045 %       mpopt = mpoption(mpopt, 'opf.dc.solver', 'CPLEX', 'verbose', 2);
0046 %
0047 %The currently defined options are as follows:
0048 %
0049 %   name                    default     description [options]
0050 %----------------------    ---------   ----------------------------------
0051 %Model options:
0052 %   model                   'AC'        AC vs. DC power flow model
0053 %       [ 'AC' - use nonlinear AC model & corresponding algorithms/options  ]
0054 %       [ 'DC' - use linear DC model & corresponding algorithms/options     ]
0055 %
0056 %Power Flow options:
0057 %   pf.alg                  'NR'        AC power flow algorithm
0058 %       [ 'NR'    - Newton's method (formulation depends on values of       ]
0059 %       [           pf.current_balance and pf.v_cartesian options)          ]
0060 %       [ 'NR-SP' - Newton's method (power mismatch, polar)                 ]
0061 %       [ 'NR-SC' - Newton's method (power mismatch, cartesian)             ]
0062 %       [ 'NR-SH' - Newton's method (power mismatch, hybrid)                ]
0063 %       [ 'NR-IP' - Newton's method (current mismatch, polar)               ]
0064 %       [ 'NR-IC' - Newton's method (current mismatch, cartesian)           ]
0065 %       [ 'NR-IH' - Newton's method (current mismatch, hybrid)              ]
0066 %       [ 'FDXB'  - Fast-Decoupled (XB version)                             ]
0067 %       [ 'FDBX'  - Fast-Decoupled (BX version)                             ]
0068 %       [ 'GS'    - Gauss-Seidel                                            ]
0069 %       [ 'PQSUM' - Power Summation method (radial networks only)           ]
0070 %       [ 'ISUM'  - Current Summation method (radial networks only)         ]
0071 %       [ 'YSUM'  - Admittance Summation method (radial networks only)      ]
0072 %   pf.current_balance     0           type of nodal balance equation
0073 %       [ 0 - use complex power balance equations                           ]
0074 %       [ 1 - use complex current balance equations                         ]
0075 %   pf.v_cartesian         0           voltage representation
0076 %       [ 0 - bus voltage variables represented in polar coordinates        ]
0077 %       [ 1 - bus voltage variables represented in cartesian coordinates    ]
0078 %       [ 2 - hybrid, polar updates computed via modified cartesian Jacobian]
0079 %   pf.tol                  1e-8        termination tolerance on per unit
0080 %                                       P & Q mismatch
0081 %   pf.nr.max_it            10          maximum number of iterations for
0082 %                                       Newton's method
0083 %   pf.nr.lin_solver        ''          linear solver passed to MPLINSOLVE to
0084 %                                       solve Newton update step
0085 %       [ ''      - default to '\' for small systems, 'LU3' for larger ones ]
0086 %       [ '\'     - built-in backslash operator                             ]
0087 %       [ 'LU'    - explicit default LU decomposition and back substitution ]
0088 %       [ 'LU3'   - 3 output arg form of LU, Gilbert-Peierls algorithm with ]
0089 %       [           approximate minimum degree (AMD) reordering             ]
0090 %       [ 'LU4'   - 4 output arg form of LU, UMFPACK solver (same as 'LU')  ]
0091 %       [ 'LU5'   - 5 output arg form of LU, UMFPACK solver, w/row scaling  ]
0092 %       [       (see MPLINSOLVE for complete list of all options)           ]
0093 %   pf.fd.max_it            30          maximum number of iterations for
0094 %                                       fast decoupled method
0095 %   pf.gs.max_it            1000        maximum number of iterations for
0096 %                                       Gauss-Seidel method
0097 %   pf.radial.max_it        20          maximum number of iterations for
0098 %                                       radial power flow methods
0099 %   pf.radial.vcorr         0           perform voltage correction procedure
0100 %                                       in distribution power flow
0101 %       [  0 - do NOT perform voltage correction                            ]
0102 %       [  1 - perform voltage correction                                   ]
0103 %   pf.enforce_q_lims       0           enforce gen reactive power limits at
0104 %                                       expense of |V|
0105 %       [  0 - do NOT enforce limits                                        ]
0106 %       [  1 - enforce limits, simultaneous bus type conversion             ]
0107 %       [  2 - enforce limits, one-at-a-time bus type conversion            ]
0108 %
0109 %Continuation Power Flow options:
0110 %   cpf.parameterization    3           choice of parameterization
0111 %       [  1 - natural                                                      ]
0112 %       [  2 - arc length                                                   ]
0113 %       [  3 - pseudo arc length                                            ]
0114 %   cpf.stop_at             'NOSE'      determins stopping criterion
0115 %       [ 'NOSE'     - stop when nose point is reached                      ]
0116 %       [ 'FULL'     - trace full nose curve                                ]
0117 %       [ <lam_stop> - stop upon reaching specified target lambda value     ]
0118 %   cpf.enforce_p_lims      0           enforce gen active power limits
0119 %       [  0 - do NOT enforce limits                                        ]
0120 %       [  1 - enforce limits, simultaneous bus type conversion             ]
0121 %   cpf.enforce_q_lims      0           enforce gen reactive power limits at
0122 %                                       expense of |V|
0123 %       [  0 - do NOT enforce limits                                        ]
0124 %       [  1 - enforce limits, simultaneous bus type conversion             ]
0125 %   cpf.enforce_v_lims      0           enforce bus voltage magnitude limits
0126 %       [  0 - do NOT enforce limits                                        ]
0127 %       [  1 - enforce limits, termination on detection                     ]
0128 %   cpf.enforce_flow_lims   0           enforce branch flow MVA limits
0129 %       [  0 - do NOT enforce limits                                        ]
0130 %       [  1 - enforce limits, termination on detection                     ]
0131 %   cpf.step                0.05        continuation power flow step size
0132 %   cpf.adapt_step          0           toggle adaptive step size feature
0133 %       [  0 - adaptive step size disabled                                  ]
0134 %       [  1 - adaptive step size enabled                                   ]
0135 %   cpf.step_min            1e-4        minimum allowed step size
0136 %   cpf.step_max            0.2         maximum allowed step size
0137 %   cpf.adapt_step_damping  0.7         damping factor for adaptive step
0138 %                                       sizing
0139 %   cpf.adapt_step_tol      1e-3        tolerance for adaptive step sizing
0140 %   cpf.target_lam_tol      1e-5        tolerance for target lambda detection
0141 %   cpf.nose_tol            1e-5        tolerance for nose point detection (pu)
0142 %   cpf.p_lims_tol          0.01        tolerance for generator active
0143 %                                       power limit enforcement (MW)
0144 %   cpf.q_lims_tol          0.01        tolerance for generator reactive
0145 %                                       power limit enforcement (MVAR)
0146 %   cpf.v_lims_tol          1e-4        tolerance for bus voltage
0147 %                                       magnitude enforcement (p.u)
0148 %   cpf.flow_lims_tol       0.01        tolerance for line MVA flow
0149 %                                       enforcement (MVA)
0150 %   cpf.plot.level          0           control plotting of noze curve
0151 %       [  0 - do not plot nose curve                                       ]
0152 %       [  1 - plot when completed                                          ]
0153 %       [  2 - plot incrementally at each iteration                         ]
0154 %       [  3 - same as 2, with 'pause' at each iteration                    ]
0155 %   cpf.plot.bus            <empty>     index of bus whose voltage is to be
0156 %                                       plotted
0157 %   cpf.user_callback       <empty>     string containing the name of a user
0158 %                                       callback function, or struct with
0159 %                                       function name, and optional priority
0160 %                                       and/or args, or cell array of such
0161 %                                       strings and/or structs, see
0162 %                                       'help cpf_default_callback' for details
0163 %
0164 %Optimal Power Flow options:
0165 %   name                    default     description [options]
0166 %----------------------    ---------   ----------------------------------
0167 %   opf.ac.solver           'DEFAULT'   AC optimal power flow solver
0168 %       [ 'DEFAULT' - choose default solver, i.e. 'MIPS'                    ]
0169 %       [ 'MIPS'    - MIPS, MATPOWER Interior Point Solver, primal/dual     ]
0170 %       [             interior point method (pure MATLAB/Octave)            ]
0171 %       [ 'FMINCON' - MATLAB Optimization Toolbox, FMINCON                  ]
0172 %       [ 'IPOPT'   - IPOPT, requires MEX interface to IPOPT solver         ]
0173 %       [             available from:                                       ]
0174 %       [                   https://github.com/coin-or/Ipopt                ]
0175 %       [ 'KNITRO'  - Artelys Knitro, requires Artelys Knitro solver,       ]
0176 %       [             available from:https://www.artelys.com/solvers/knitro/]
0177 %       [ 'MINOPF'  - MINOPF, MINOS-based solver, requires optional         ]
0178 %       [             MEX-based MINOPF package, available from:             ]
0179 %       [                   http://www.pserc.cornell.edu/minopf/            ]
0180 %       [ 'PDIPM'   - PDIPM, primal/dual interior point method, requires    ]
0181 %       [             optional MEX-based TSPOPF package, available from:    ]
0182 %       [                   http://www.pserc.cornell.edu/tspopf/            ]
0183 %       [ 'SDPOPF'  - SDPOPF, solver based on semidefinite relaxation of    ]
0184 %       [             OPF problem, requires optional packages:              ]
0185 %       [               SDP_PF, available in extras/sdp_pf                  ]
0186 %       [               YALMIP, available from:                             ]
0187 %       [                   https://yalmip.github.io                        ]
0188 %       [               SDP solver such as SeDuMi, available from:          ]
0189 %       [                   http://sedumi.ie.lehigh.edu/                    ]
0190 %       [ 'TRALM'   - TRALM, trust region based augmented Langrangian       ]
0191 %       [             method, requires TSPOPF (see 'PDIPM')                 ]
0192 %   opf.dc.solver           'DEFAULT'   DC optimal power flow solver
0193 %       [ 'DEFAULT' - choose solver based on availability in the following  ]
0194 %       [             order: 'GUROBI', 'CPLEX', 'MOSEK', 'OT',              ]
0195 %       [             'GLPK' (linear costs only), 'BPMPD', 'MIPS'           ]
0196 %       [ 'MIPS'    - MIPS, MATPOWER Interior Point Solver, primal/dual     ]
0197 %       [             interior point method (pure MATLAB/Octave)            ]
0198 %       [ 'BPMPD'   - BPMPD, requires optional MEX-based BPMPD_MEX package  ]
0199 %       [             available from: http://www.pserc.cornell.edu/bpmpd/   ]
0200 %       [ 'CLP'     - CLP, requires interface to COIN-OP LP solver          ]
0201 %       [             available from:https://github.com/coin-or/Clp         ]
0202 %       [ 'CPLEX'   - CPLEX, requires CPLEX solver available from:          ]
0203 %       [             https://www.ibm.com/analytics/cplex-optimizer         ]
0204 %       [ 'GLPK'    - GLPK, requires interface to GLPK solver               ]
0205 %       [             available from: https://www.gnu.org/software/glpk/    ]
0206 %       [             (GLPK does not work with quadratic cost functions)    ]
0207 %       [ 'GUROBI'  - GUROBI, requires Gurobi optimizer (v. 5+)             ]
0208 %       [             available from: https://www.gurobi.com/               ]
0209 %       [ 'IPOPT'   - IPOPT, requires MEX interface to IPOPT solver         ]
0210 %       [             available from:                                       ]
0211 %       [                 https://github.com/coin-or/Ipopt                  ]
0212 %       [ 'MOSEK'   - MOSEK, requires MATLAB interface to MOSEK solver      ]
0213 %       [             available from: https://www.mosek.com/                ]
0214 %       [ 'OSQP'    - OSQP, requires MATLAB interface to OSQP solver        ]
0215 %       [             available from: https://osqp.org/                     ]
0216 %       [ 'OT'      - MATLAB Optimization Toolbox, QUADPROG, LINPROG        ]
0217 %   opf.current_balance     0           type of nodal balance equation
0218 %       [ 0 - use complex power balance equations                           ]
0219 %       [ 1 - use complex current balance equations                         ]
0220 %   opf.v_cartesian         0           voltage representation
0221 %       [ 0 - bus voltage variables represented in polar coordinates        ]
0222 %       [ 1 - bus voltage variables represented in cartesian coordinates    ]
0223 %   opf.violation           5e-6        constraint violation tolerance
0224 %   opf.use_vg              0           respect gen voltage setpt     [ 0-1 ]
0225 %       [ 0 - use specified bus Vmin & Vmax, and ignore gen Vg              ]
0226 %       [ 1 - replace specified bus Vmin & Vmax by corresponding gen Vg     ]
0227 %       [ between 0 and 1 - use a weighted average of the 2 options         ]
0228 %   opf.flow_lim            'S'         quantity limited by branch flow
0229 %                                       constraints
0230 %       [ 'S' - apparent power flow (limit in MVA)                          ]
0231 %       [ 'P' - active power flow, implemented using P   (limit in MW)      ]
0232 %       [ '2' - active power flow, implemented using P^2 (limit in MW)      ]
0233 %       [ 'I' - current magnitude (limit in MVA at 1 p.u. voltage)          ]
0234 %   opf.ignore_angle_lim    0           angle diff limits for branches
0235 %       [ 0 - include angle difference limits, if specified                 ]
0236 %       [ 1 - ignore angle difference limits even if specified              ]
0237 %   opf.softlims.default    1           behavior of OPF soft limits for
0238 %                                       which parameters are not explicitly
0239 %                                       provided
0240 %       [ 0 - do not include softlims if not explicitly specified           ]
0241 %       [ 1 - include softlims w/default values if not explicitly specified ]
0242 %   opf.init_from_mpc       -1          (DEPRECATED: use opf.start instead)
0243 %                                       specify whether to use current state
0244 %                                       in MATPOWER case to initialize OPF
0245 %                                       (currently only supported by fmincon,
0246 %                                        Ipopt, Knitro and MIPS solvers,
0247 %                                        others always use mpc)
0248 %       [  -1 - MATPOWER decides, based on solver/algorithm                 ]
0249 %       [   0 - ignore current state in MATPOWER case (only applies to      ]
0250 %       [       fmincon, Ipopt, Knitro and MIPS, which use an interior pt   ]
0251 %       [       estimate; others use current state as with opf.start = 2)   ]
0252 %       [   1 - use current state in MATPOWER case                          ]
0253 %   opf.start               0           strategy for initializing OPF start pt
0254 %       [   0 - default, MATPOWER decides based on solver                   ]
0255 %       [       (currently identical to 1)                                  ]
0256 %       [   1 - ignore current state in MATPOWER case (only applies to      ]
0257 %       [       fmincon, Ipopt, Knitro and MIPS, which use an interior pt   ]
0258 %       [       estimate; others use current state as with opf.start = 2)   ]
0259 %       [   2 - use current state in MATPOWER case                          ]
0260 %       [   3 - solve power flow and use resulting state                    ]
0261 %   opf.return_raw_der      0           for AC OPF, return constraint and
0262 %                                       derivative info in results.raw
0263 %                                       (in fields g, dg, df, d2f) [ 0 or 1 ]
0264 %
0265 %Output options:
0266 %   name                    default     description [options]
0267 %----------------------    ---------   ----------------------------------
0268 %   verbose                 1           amount of progress info printed
0269 %       [   0 - print no progress info                                      ]
0270 %       [   1 - print a little progress info                                ]
0271 %       [   2 - print a lot of progress info                                ]
0272 %       [   3 - print all progress info                                     ]
0273 %   out.all                 -1          controls pretty-printing of results
0274 %       [  -1 - individual flags control what prints                        ]
0275 %       [   0 - do not print anything (overrides individual flags, ignored  ]
0276 %       [       for files specified as FNAME arg to runpf(), runopf(), etc.)]
0277 %       [   1 - print everything (overrides individual flags)               ]
0278 %   out.sys_sum             1           print system summary       [ 0 or 1 ]
0279 %   out.area_sum            0           print area summaries       [ 0 or 1 ]
0280 %   out.bus                 1           print bus detail           [ 0 or 1 ]
0281 %   out.branch              1           print branch detail        [ 0 or 1 ]
0282 %   out.gen                 0           print generator detail     [ 0 or 1 ]
0283 %   out.lim.all             -1          controls constraint info output
0284 %       [  -1 - individual flags control what constraint info prints        ]
0285 %       [   0 - no constraint info (overrides individual flags)             ]
0286 %       [   1 - binding constraint info (overrides individual flags)        ]
0287 %       [   2 - all constraint info (overrides individual flags)            ]
0288 %   out.lim.v               1           control voltage limit info
0289 %       [   0 - do not print                                                ]
0290 %       [   1 - print binding constraints only                              ]
0291 %       [   2 - print all constraints                                       ]
0292 %       [   (same options for OUT_LINE_LIM, OUT_PG_LIM, OUT_QG_LIM)         ]
0293 %   out.lim.line            1           control line flow limit info
0294 %   out.lim.pg              1           control gen active power limit info
0295 %   out.lim.qg              1           control gen reactive pwr limit info
0296 %   out.force               0           print results even if success
0297 %                                       flag = 0                   [ 0 or 1 ]
0298 %   out.suppress_detail     -1          suppress all output but system summary
0299 %       [  -1 - suppress details for large systems (> 500 buses)            ]
0300 %       [   0 - do not suppress any output specified by other flags         ]
0301 %       [   1 - suppress all output except system summary section           ]
0302 %       [       (overrides individual flags, but not out.all = 1)           ]
0303 %
0304 %Solver specific options:
0305 %       name                    default     description [options]
0306 %   -----------------------    ---------   ----------------------------------
0307 %   MIPS:
0308 %       mips.linsolver          ''          linear system solver
0309 %           [   '' or '\'   build-in backslash \ operator (e.g. x = A \ b)  ]
0310 %           [   'PARDISO'   PARDISO solver (if available)                   ]
0311 %       mips.feastol            0           feasibility (equality) tolerance
0312 %                                           (set to opf.violation by default)
0313 %       mips.gradtol            1e-6        gradient tolerance
0314 %       mips.comptol            1e-6        complementary condition
0315 %                                           (inequality) tolerance
0316 %       mips.costtol            1e-6        optimality tolerance
0317 %       mips.max_it             150         maximum number of iterations
0318 %       mips.step_control       0           enable step-size cntrl [ 0 or 1 ]
0319 %       mips.sc.red_it          20          maximum number of reductions per
0320 %                                           iteration with step control
0321 %       mips.xi                 0.99995     constant used in alpha updates*
0322 %       mips.sigma              0.1         centering parameter*
0323 %       mips.z0                 1           used to initialize slack variables*
0324 %       mips.alpha_min          1e-8        returns "Numerically Failed" if
0325 %                                           either alpha parameter becomes
0326 %                                           smaller than this value*
0327 %       mips.rho_min            0.95        lower bound on rho_t*
0328 %       mips.rho_max            1.05        upper bound on rho_t*
0329 %       mips.mu_threshold       1e-5        KT multipliers smaller than this
0330 %                                           value for non-binding constraints
0331 %                                           are forced to zero
0332 %       mips.max_stepsize       1e10        returns "Numerically Failed" if the
0333 %                                           2-norm of the reduced Newton step
0334 %                                           exceeds this value*
0335 %           * See the corresponding Appendix in the manual for details.
0336 %
0337 %   CPLEX:
0338 %       cplex.lpmethod          0           solution algorithm for LP problems
0339 %           [   0 - automatic: let CPLEX choose                             ]
0340 %           [   1 - primal simplex                                          ]
0341 %           [   2 - dual simplex                                            ]
0342 %           [   3 - network simplex                                         ]
0343 %           [   4 - barrier                                                 ]
0344 %           [   5 - sifting                                                 ]
0345 %           [   6 - concurrent (dual, barrier, and primal)                  ]
0346 %       cplex.qpmethod          0           solution algorithm for QP problems
0347 %           [   0 - automatic: let CPLEX choose                             ]
0348 %           [   1 - primal simplex optimizer                                ]
0349 %           [   2 - dual simplex optimizer                                  ]
0350 %           [   3 - network optimizer                                       ]
0351 %           [   4 - barrier optimizer                                       ]
0352 %       cplex.opts              <empty>     see CPLEX_OPTIONS for details
0353 %       cplex.opt_fname         <empty>     see CPLEX_OPTIONS for details
0354 %       cplex.opt               0           see CPLEX_OPTIONS for details
0355 %
0356 %   FMINCON:
0357 %       fmincon.alg             4           algorithm used by fmincon() for OPF
0358 %                                           for Opt Toolbox 4 and later
0359 %            [  1 - active-set (not suitable for large problems)            ]
0360 %            [  2 - interior-point, w/default 'bfgs' Hessian approx         ]
0361 %            [  3 - interior-point, w/ 'lbfgs' Hessian approx               ]
0362 %            [  4 - interior-point, w/exact user-supplied Hessian           ]
0363 %            [  5 - interior-point, w/Hessian via finite differences        ]
0364 %            [  6 - sqp (not suitable for large problems)                   ]
0365 %       fmincon.tol_x           1e-4        termination tol on x
0366 %       fmincon.tol_f           1e-4        termination tol on f
0367 %       fmincon.max_it          0           maximum number of iterations
0368 %                                                           [  0 => default ]
0369 %
0370 %   GUROBI:
0371 %       gurobi.method           0           solution algorithm (Method)
0372 %           [  -1 - automatic, let Gurobi decide                            ]
0373 %           [   0 - primal simplex                                          ]
0374 %           [   1 - dual simplex                                            ]
0375 %           [   2 - barrier                                                 ]
0376 %           [   3 - concurrent (LP only)                                    ]
0377 %           [   4 - deterministic concurrent (LP only)                      ]
0378 %       gurobi.timelimit        Inf         maximum time allowed (TimeLimit)
0379 %       gurobi.threads          0           max number of threads (Threads)
0380 %       gurobi.opts             <empty>     see GUROBI_OPTIONS for details
0381 %       gurobi.opt_fname        <empty>     see GUROBI_OPTIONS for details
0382 %       gurobi.opt              0           see GUROBI_OPTIONS for details
0383 %
0384 %   IPOPT:
0385 %       ipopt.opts              <empty>     see IPOPT_OPTIONS for details
0386 %       ipopt.opt_fname         <empty>     see IPOPT_OPTIONS for details
0387 %       ipopt.opt               0           see IPOPT_OPTIONS for details
0388 %
0389 %   KNITRO:
0390 %       knitro.tol_x            1e-4        termination tol on x
0391 %       knitro.tol_f            1e-4        termination tol on f
0392 %       knitro.maxit            0           maximum number of iterations
0393 %                                                           [  0 => default ]
0394 %       knitro.opt_fname        <empty>     name of user-supplied native
0395 %                                           Knitro options file that overrides
0396 %                                           all other options
0397 %       knitro.opt              0           if knitro.opt_fname is empty and
0398 %                                           knitro.opt is a non-zero integer N
0399 %                                           then knitro.opt_fname is auto-
0400 %                                           generated as:
0401 %                                           'knitro_user_options_N.txt'
0402 %
0403 %   LINPROG:
0404 %       linprog                 <empty>     LINPROG options passed to
0405 %                                           OPTIMOPTIONS or OPTIMSET.
0406 %                                           see LINPROG in the Optimization
0407 %                                           Toolbox for details
0408 %
0409 %   MINOPF:
0410 %       minopf.feastol          0 (1e-3)    primal feasibility tolerance
0411 %                                           (set to opf.violation by default)
0412 %       minopf.rowtol           0 (1e-3)    row tolerance
0413 %       minopf.xtol             0 (1e-4)    x tolerance
0414 %       minopf.majdamp          0 (0.5)     major damping parameter
0415 %       minopf.mindamp          0 (2.0)     minor damping parameter
0416 %       minopf.penalty          0 (1.0)     penalty parameter
0417 %       minopf.major_it         0 (200)     major iterations
0418 %       minopf.minor_it         0 (2500)    minor iterations
0419 %       minopf.max_it           0 (2500)    iterations limit
0420 %       minopf.verbosity        -1          amount of progress info printed
0421 %           [  -1 - controlled by 'verbose' option                          ]
0422 %           [   0 - print nothing                                           ]
0423 %           [   1 - print only termination status message                   ]
0424 %           [   2 - print termination status and screen progress            ]
0425 %           [   3 - print screen progress, report file (usually fort.9)     ]
0426 %       minopf.core             0 (1200*nb + 2*(nb+ng)^2) memory allocation
0427 %       minopf.supbasic_lim     0 (2*nb + 2*ng) superbasics limit
0428 %       minopf.mult_price       0 (30)      multiple price
0429 %
0430 %   MOSEK:
0431 %       mosek.lp_alg            0           solution algorithm
0432 %                                               (MSK_IPAR_OPTIMIZER)
0433 %           for MOSEK 8.x ...        (see MOSEK_SYMBCON for a "better way")
0434 %           [   0 - automatic: let MOSEK choose                             ]
0435 %           [   1 - dual simplex                                            ]
0436 %           [   2 - automatic: let MOSEK choose                             ]
0437 %           [   3 - automatic simplex (MOSEK chooses which simplex method)  ]
0438 %           [   4 - interior point                                          ]
0439 %           [   6 - primal simplex                                          ]
0440 %       mosek.max_it            0 (400)     interior point max iterations
0441 %                                               (MSK_IPAR_INTPNT_MAX_ITERATIONS)
0442 %       mosek.gap_tol           0 (1e-8)    interior point relative gap tol
0443 %                                               (MSK_DPAR_INTPNT_TOL_REL_GAP)
0444 %       mosek.max_time          0 (-1)      maximum time allowed
0445 %                                               (MSK_DPAR_OPTIMIZER_MAX_TIME)
0446 %       mosek.num_threads       0 (1)       max number of threads
0447 %                                               (MSK_IPAR_INTPNT_NUM_THREADS)
0448 %       mosek.opts              <empty>     see MOSEK_OPTIONS for details
0449 %       mosek.opt_fname         <empty>     see MOSEK_OPTIONS for details
0450 %       mosek.opt               0           see MOSEK_OPTIONS for details
0451 %
0452 %   OSQP:
0453 %       osqp.opts               <empty>     see OSQP_OPTIONS for details
0454 %
0455 %   QUADPROG:
0456 %       quadprog                <empty>     QUADPROG options passed to
0457 %                                           OPTIMOPTIONS or OPTIMSET.
0458 %                                           see QUADPROG in the Optimization
0459 %                                           Toolbox for details
0460 %
0461 %   TSPOPF:
0462 %       pdipm.feastol           0           feasibility (equality) tolerance
0463 %                                           (set to opf.violation by default)
0464 %       pdipm.gradtol           1e-6        gradient tolerance
0465 %       pdipm.comptol           1e-6        complementary condition
0466 %                                           (inequality) tolerance
0467 %       pdipm.costtol           1e-6        optimality tolerance
0468 %       pdipm.max_it            150         maximum number of iterations
0469 %       pdipm.step_control      0           enable step-size cntrl [ 0 or 1 ]
0470 %       pdipm.sc.red_it         20          maximum number of reductions per
0471 %                                           iteration with step control
0472 %       pdipm.sc.smooth_ratio   0.04        piecewise linear curve smoothing
0473 %                                           ratio
0474 %
0475 %       tralm.feastol           0           feasibility tolerance
0476 %                                           (set to opf.violation by default)
0477 %       tralm.primaltol         5e-4        primal variable tolerance
0478 %       tralm.dualtol           5e-4        dual variable tolerance
0479 %       tralm.costtol           1e-5        optimality tolerance
0480 %       tralm.major_it          40          maximum number of major iterations
0481 %       tralm.minor_it          40          maximum number of minor iterations
0482 %       tralm.smooth_ratio      0.04        piecewise linear curve smoothing
0483 %                                           ratio
0484 %
0485 %Experimental Options:
0486 %   exp.sys_wide_zip_loads.pw   <empty>     1 x 3 vector of active load fraction
0487 %                                           to be modeled as constant power,
0488 %                                           constant current and constant
0489 %                                           impedance, respectively, where
0490 %                                           <empty> means use [1 0 0]
0491 %   exp.sys_wide_zip_loads.qw   <empty>     same for reactive power, where
0492 %                                           <empty> means use same value as
0493 %                                           for 'pw'
0494 
0495 %   MATPOWER
0496 %   Copyright (c) 2013-2017, Power Systems Engineering Research Center (PSERC)
0497 %   by Ray Zimmerman, PSERC Cornell
0498 %   and Shrirang Abhyankar, Argonne National Laboratory
0499 %
0500 %   This file is part of MATPOWER.
0501 %   Covered by the 3-clause BSD License (see LICENSE file for details).
0502 %   See https://matpower.org for more info.
0503 
0504 %% some constants
0505 N = 124;                %% number of options in old-style vector (MATPOWER 4.1)
0506 N40 = 116;              %% dimension of MATPOWER 4.0 options vector
0507 N32 = 93;               %% dimension of MATPOWER 3.2 options vector
0508 v = mpoption_version;   %% version number of MATPOWER options struct
0509 
0510 %% cache nested_struct_copy options to speed things up
0511 persistent nsc_opt;
0512 
0513 %% initialize flags and arg counter
0514 have_opt0 = 0;          %% existing options struct or vector provided?
0515 have_old_style_ov = 0;  %% override options using old-style names?
0516 return_old_style = 0;   %% return value as old-style vector?
0517 k = 1;
0518 if nargin > 0
0519     opt0 = varargin{k};
0520     if isstruct(opt0) && isfield(opt0, 'v')         %% options struct
0521         have_opt0 = 1;
0522         k = k + 1;
0523     elseif isnumeric(opt0) && size(opt0, 2) == 1    %% options vector
0524         nn = size(opt0, 1);
0525         if ismember(nn, [N N40 N32])                %% of valid size
0526             %% expand smaller option vectors (from before MATPOWER 4.1)
0527             if nn < N
0528                 optv = mpoption_old();
0529                 opt0(nn+1:N) = optv(nn+1:N);
0530             end
0531             have_opt0 = 1;
0532             k = k + 1;
0533         end
0534     end
0535 end
0536 
0537 %% create base options vector to which overrides are made
0538 if have_opt0
0539     if isstruct(opt0)               %% it's already a valid options struct
0540         if DEBUG, fprintf('OPT0 is a valid options struct\n'); end
0541         if opt0.v < v
0542             %% convert older version to current version
0543             opt_d = mpoption_default();
0544             if opt0.v == 1          %% convert version 1 to 2
0545                 if isfield(opt_d, 'linprog')
0546                     opt0.lingprog = opt_d.linprog;
0547                 end
0548                 if isfield(opt_d, 'quadprog')
0549                     opt0.quadprog = opt_d.quadprog;
0550                 end
0551             end
0552             if opt0.v <= 2          %% convert version 2 to 3
0553                 opt0.out.suppress_detail = opt_d.out.suppress_detail;
0554             end
0555             %if opt0.v <= 3          %% convert version 3 to 4
0556                 %% new mips options were all optional, no conversion needed
0557             %end
0558             if opt0.v <= 4          %% convert version 4 to 5
0559                 opt0.opf.init_from_mpc = opt_d.opf.init_from_mpc;
0560             end
0561             if opt0.v <= 5          %% convert version 5 to 6
0562                 if isfield(opt_d, 'clp')
0563                     opt0.clp = opt_d.clp;
0564                 end
0565             end
0566             if opt0.v <= 6          %% convert version 6 to 7
0567                 if isfield(opt_d, 'intlinprog')
0568                     opt0.intlinprog = opt_d.intlinprog;
0569                 end
0570             end
0571             if opt0.v <= 7          %% convert version 7 to 8
0572                 opt0.mips.linsolver = opt_d.mips.linsolver;
0573             end
0574             if opt0.v <= 8          %% convert version 8 to 9
0575                 opt0.exp.sys_wide_zip_loads = opt_d.exp.sys_wide_zip_loads;
0576             end
0577             if opt0.v <= 9          %% convert version 9 to 10
0578                 opt0.most = opt_d.most;
0579             end
0580             if opt0.v <= 10         %% convert version 10 to 11
0581                 opt0.cpf.enforce_p_lims     = opt_d.cpf.enforce_p_lims;
0582                 opt0.cpf.enforce_q_lims     = opt_d.cpf.enforce_q_lims;
0583                 opt0.cpf.adapt_step_damping = opt_d.cpf.adapt_step_damping;
0584                 opt0.cpf.target_lam_tol     = opt_d.cpf.target_lam_tol;
0585                 opt0.cpf.nose_tol           = opt_d.cpf.nose_tol;
0586                 opt0.cpf.p_lims_tol         = opt_d.cpf.p_lims_tol;
0587                 opt0.cpf.q_lims_tol         = opt_d.cpf.q_lims_tol;
0588                 if (~isempty(opt0.cpf.user_callback_args) && ...
0589                         ~isstruct(opt0.cpf.user_callback_args)) || ...
0590                         (isstruct(opt0.cpf.user_callback_args) && ...
0591                          ~isempty(fields(opt0.cpf.user_callback_args)))
0592                     warning('The ''cpf.user_callback_args'' option has been removed. Please include the args in a struct in ''cpf.user_callback'' instead.')
0593                 end
0594                 opt0.cpf = rmfield(opt0.cpf, 'user_callback_args');
0595             end
0596             if opt0.v <= 11         %% convert version 11 to 12
0597                 opt0.opf.use_vg = opt_d.opf.use_vg;
0598             end
0599             if opt0.v <= 12         %% convert version 12 to 13
0600                 opt0.pf.radial.max_it   = opt_d.pf.radial.max_it;
0601                 opt0.pf.radial.vcorr    = opt_d.pf.radial.vcorr;
0602             end
0603             if opt0.v <= 13         %% convert version 13 to 14
0604                 opt0.pf.nr.lin_solver   = opt_d.pf.nr.lin_solver;
0605             end
0606             if opt0.v <= 14         %% convert version 14 to 15
0607                 opt0.cpf.enforce_v_lims     = opt_d.cpf.enforce_v_lims;
0608                 opt0.cpf.enforce_flow_lims  = opt_d.cpf.enforce_flow_lims;
0609                 opt0.cpf.v_lims_tol         = opt_d.cpf.v_lims_tol;
0610                 opt0.cpf.flow_lims_tol      = opt_d.cpf.flow_lims_tol;
0611             end
0612             if opt0.v <= 15          %% convert version 15 to 16
0613                 opt0.opf.start  = opt_d.opf.start;
0614             end
0615             if opt0.v <= 16          %% convert version 16 to 17
0616                 if isfield(opt_d, 'knitro')
0617                     opt0.knitro.maxit  = opt_d.knitro.maxit;
0618                 end
0619             end
0620             if opt0.v <= 17         %% convert version 17 to 18
0621                 opt0.opf.current_balance = opt_d.opf.current_balance;
0622                 opt0.opf.v_cartesian     = opt_d.opf.v_cartesian;
0623             end
0624             if opt0.v <= 18         %% convert version 18 to 19
0625                 opt0.opf.softlims.default = opt_d.opf.softlims.default;
0626             end
0627             if opt0.v <= 19         %% convert version 19 to 20
0628                 opt0.pf.current_balance = opt_d.pf.current_balance;
0629                 opt0.pf.v_cartesian     = opt_d.pf.v_cartesian;
0630             end
0631             if opt0.v <= 20         %% convert version 20 to 21
0632                 opt0.osqp = opt_d.osqp;
0633             end
0634             opt0.v = v;
0635         end
0636         opt = opt0;
0637     else                            %% convert from old-style options vector
0638         if DEBUG, fprintf('OPT0 is a old-style options vector\n'); end
0639         opt = mpoption_v2s(opt0);
0640     end
0641 else                                %% use default options struct as base
0642     if DEBUG, fprintf('no OPT0, starting with default options struct\n'); end
0643     opt = mpoption_default();
0644 end
0645 
0646 
0647 %% do we have OVERRIDES or NAME/VALUE pairs
0648 ov = [];
0649 if nargin - k == 0          %% looking at last arg, must be OVERRIDES
0650     if isstruct(varargin{k})        %% OVERRIDES provided as struct
0651         if DEBUG, fprintf('OVERRIDES struct\n'); end
0652         ov = varargin{k};
0653     elseif ischar(varargin{k})      %% OVERRIDES provided as file/function name
0654         if DEBUG, fprintf('OVERRIDES file/function name\n'); end
0655         try
0656             ov = feval(varargin{k});
0657         catch
0658             error('mpoption: Unable to load MATPOWER options from ''%s''', varargin{k});
0659         end
0660         if ~isstruct(ov)
0661             error('mpoption: calling ''%s'' did not return a struct', varargin{k});
0662         end
0663     elseif isempty(varargin{k})
0664         return_old_style = 1;
0665     else
0666         error('mpoption: OVERRIDES must be a struct or the name of a function that returns a struct');
0667     end
0668 elseif nargin - k > 0 && mod(nargin-k, 2)   %% even number of remaining args
0669     if DEBUG, fprintf('NAME/VALUE pairs override defaults\n'); end
0670     %% process NAME/VALUE pairs
0671     if strcmp(varargin{k}, upper(varargin{k}))  %% old-style, all UPPERCASE option pairs
0672         %% NOTE: new top-level option fields cannot be all uppercase
0673         if ~have_opt0
0674             opt_v = mpoption_old(varargin{:});  %% create modified vector ...
0675             opt = mpoption_v2s(opt_v);          %% ... then convert
0676         else
0677             have_old_style_ov = 1;
0678             %% convert pairs to struct
0679             while k < nargin
0680                 name = varargin{k};
0681                 val  = varargin{k+1};
0682                 k = k + 2;
0683                 ov.(name) = val;
0684             end
0685         end
0686     else                                        %% new option pairs
0687         %% convert pairs to struct
0688         while k < nargin
0689             name = varargin{k};
0690             val  = varargin{k+1};
0691             k = k + 2;
0692             c = regexp(name, '([^\.]*)', 'tokens');
0693             s = struct();
0694             for i = 1:length(c)
0695                 s(i).type = '.';
0696                 s(i).subs = c{i}{1};
0697             end
0698             ov = subsasgn(ov, s, val);
0699         end
0700     end
0701 elseif nargin == 0 || nargin == 1
0702     if DEBUG, fprintf('no OVERRIDES, return default options struct or converted OPT0 vector\n'); end
0703 else
0704     error('mpoption: invalid calling syntax, see ''help mpoption'' to double-check the valid options');
0705 end
0706 
0707 %% apply overrides
0708 if ~isempty(ov)
0709     if have_old_style_ov
0710         opt = apply_old_mpoption_overrides(opt, ov);
0711     else
0712         if ~isstruct(nsc_opt)   %% nsc_opt is cached to speed things up
0713             vf = nested_struct_copy(mpoption_default(), mpoption_info_mips('V'));
0714             vf = nested_struct_copy(vf, mpoption_optional_fields());
0715             ex = struct(...
0716                 'name', {}, ...
0717                 'check', {}, ...
0718                 'copy_mode', {} ...
0719             );
0720             %% add exceptions for optional packages
0721             opt_pkgs = mpoption_optional_pkgs();
0722             n = length(ex);
0723             for k = 1:length(opt_pkgs)
0724                 fname = ['mpoption_info_' opt_pkgs{k}];
0725                 if exist(fname, 'file') == 2
0726                     opt_ex = feval(fname, 'E');
0727                     nex = length(opt_ex);
0728                     if ~isempty(opt_ex)
0729                         for j = 1:nex
0730                             ex(n+j).name = opt_ex(j).name;
0731                         end
0732                         if isfield(opt_ex, 'check')
0733                             for j = 1:nex
0734                                 ex(n+j).check = opt_ex(j).check;
0735                             end
0736                         end
0737                         if isfield(opt_ex, 'copy_mode')
0738                             for j = 1:nex
0739                                 ex(n+j).copy_mode = opt_ex(j).copy_mode;
0740                             end
0741                         end
0742                         if isfield(opt_ex, 'valid_fields')
0743                             for j = 1:nex
0744                                 ex(n+j).valid_fields = opt_ex(j).valid_fields;
0745                             end
0746                         end
0747                         n = n + nex;
0748                     end
0749                 end
0750             end
0751             nsc_opt = struct('check', 1, 'valid_fields', vf, 'exceptions', ex);
0752         end
0753 %         if have_feature('catchme')
0754 %             try
0755 %                 opt = nested_struct_copy(opt, ov, nsc_opt);
0756 %             catch me
0757 %                 str = strrep(me.message, 'field', 'option');
0758 %                 str = strrep(str, 'nested_struct_copy', 'mpoption');
0759 %                 error(str);
0760 %             end
0761 %         else
0762             try
0763                 opt = nested_struct_copy(opt, ov, nsc_opt);
0764             catch
0765                 me = lasterr;
0766                 str = strrep(me, 'field', 'option');
0767                 str = strrep(str, 'nested_struct_copy', 'mpoption');
0768                 error(str);
0769             end
0770 %         end
0771     end
0772 end
0773 if return_old_style
0774     opt = mpoption_s2v(opt);
0775 end
0776 
0777 
0778 %%-------------------------------------------------------------------
0779 function opt = apply_old_mpoption_overrides(opt0, ov)
0780 %
0781 %   OPT0 is assumed to already have all of the fields and sub-fields found
0782 %   in the default options struct.
0783 
0784 %% initialize output
0785 opt = opt0;
0786 
0787 errstr = 'mpoption: %g is not a valid value for the old-style ''%s'' option';
0788 fields = fieldnames(ov);
0789 for f = 1:length(fields)
0790     ff = fields{f};
0791     switch ff
0792         case 'PF_ALG'
0793             switch ov.(ff)
0794                 case 1
0795                     opt.pf.alg = 'NR';      %% Newton's method
0796                 case 2
0797                     opt.pf.alg = 'FDXB';    %% fast-decoupled (XB version)
0798                 case 3
0799                     opt.pf.alg = 'FDBX';    %% fast-decoupled (BX version)
0800                 case 4
0801                     opt.pf.alg = 'GS';      %% Gauss-Seidel
0802                 otherwise
0803                     error(errstr, ov.(ff), ff);
0804             end
0805         case 'PF_TOL'
0806             opt.pf.tol = ov.(ff);
0807         case 'PF_MAX_IT'
0808             opt.pf.nr.max_it = ov.(ff);
0809         case 'PF_MAX_IT_FD'
0810             opt.pf.fd.max_it = ov.(ff);
0811         case 'PF_MAX_IT_GS'
0812             opt.pf.gs.max_it = ov.(ff);
0813         case 'ENFORCE_Q_LIMS'
0814             opt.pf.enforce_q_lims = ov.(ff);
0815         case 'PF_DC'
0816             switch ov.(ff)
0817                 case 0
0818                     opt.model = 'AC';
0819                 case 1
0820                     opt.model = 'DC';
0821                 otherwise
0822                     error(errstr, ov.(ff), ff);
0823             end
0824         case 'OPF_ALG'
0825             switch ov.(ff)
0826                 case 0
0827                     opt.opf.ac.solver = 'DEFAULT';
0828                 case 500
0829                     opt.opf.ac.solver = 'MINOPF';
0830                 case 520
0831                     opt.opf.ac.solver = 'FMINCON';
0832                 case {540, 545}
0833                     opt.opf.ac.solver = 'PDIPM';
0834                     if ov.(ff) == 545
0835                         opt.pdipm.step_control = 1;
0836                     else
0837                         opt.pdipm.step_control = 0;
0838                     end
0839                 case 550
0840                     opt.opf.ac.solver = 'TRALM';
0841                 case {560, 565}
0842                     opt.opf.ac.solver = 'MIPS';
0843                     if ov.(ff) == 565
0844                         opt.mips.step_control = 1;
0845                     else
0846                         opt.mips.step_control = 0;
0847                     end
0848                 case 580
0849                     opt.opf.ac.solver = 'IPOPT';
0850                 case 600
0851                     opt.opf.ac.solver = 'KNITRO';
0852                 otherwise
0853                     error(errstr, ov.(ff), ff);
0854             end
0855         case 'OPF_VIOLATION'
0856             opt.opf.violation = ov.(ff);
0857         case 'CONSTR_TOL_X'
0858             opt.fmincon.tol_x = ov.(ff);
0859             opt.knitro.tol_x = ov.(ff);
0860         case 'CONSTR_TOL_F'
0861             opt.fmincon.tol_f = ov.(ff);
0862             opt.knitro.tol_f = ov.(ff);
0863         case 'CONSTR_MAX_IT'
0864             opt.fmincon.max_it = ov.(ff);
0865         case 'OPF_FLOW_LIM'
0866             switch ov.(ff)
0867                 case 0
0868                     opt.opf.flow_lim = 'S';   %% apparent power (MVA)
0869                 case 1
0870                     opt.opf.flow_lim = 'P';   %% real power (MW)
0871                 case 2
0872                     opt.opf.flow_lim = 'I';   %% current magnitude (MVA @ 1 p.u. voltage)
0873                 otherwise
0874                     error(errstr, ov.(ff), ff);
0875             end
0876         case 'OPF_IGNORE_ANG_LIM'
0877             opt.opf.ignore_angle_lim = ov.(ff);
0878         case 'OPF_ALG_DC'
0879             switch ov.(ff)
0880                 case 0
0881                     opt.opf.dc.solver = 'DEFAULT';
0882                 case 100
0883                     opt.opf.dc.solver = 'BPMPD';
0884                 case {200, 250}
0885                     opt.opf.dc.solver = 'MIPS';
0886                     if ov.(ff) == 250
0887                         opt.mips.step_control = 1;
0888                     else
0889                         opt.mips.step_control = 0;
0890                     end
0891                 case 300
0892                     opt.opf.dc.solver = 'OT';     %% QUADPROG, LINPROG
0893                 case 400
0894                     opt.opf.dc.solver = 'IPOPT';
0895                 case 500
0896                     opt.opf.dc.solver = 'CPLEX';
0897                 case 600
0898                     opt.opf.dc.solver = 'MOSEK';
0899                 case 700
0900                     opt.opf.dc.solver = 'GUROBI';
0901                 otherwise
0902                     error(errstr, ov.(ff), ff);
0903             end
0904         case 'VERBOSE'
0905             opt.verbose = ov.(ff);
0906         case 'OUT_ALL'
0907             opt.out.all = ov.(ff);
0908         case 'OUT_SYS_SUM'
0909             opt.out.sys_sum = ov.(ff);
0910         case 'OUT_AREA_SUM'
0911             opt.out.area_sum = ov.(ff);
0912         case 'OUT_BUS'
0913             opt.out.bus = ov.(ff);
0914         case 'OUT_BRANCH'
0915             opt.out.branch = ov.(ff);
0916         case 'OUT_GEN'
0917             opt.out.gen = ov.(ff);
0918         case 'OUT_ALL_LIM'
0919             opt.out.lim.all = ov.(ff);
0920         case 'OUT_V_LIM'
0921             opt.out.lim.v = ov.(ff);
0922         case 'OUT_LINE_LIM'
0923             opt.out.lim.line = ov.(ff);
0924         case 'OUT_PG_LIM'
0925             opt.out.lim.pg = ov.(ff);
0926         case 'OUT_QG_LIM'
0927             opt.out.lim.qg = ov.(ff);
0928         case 'OUT_FORCE'
0929             opt.out.force = ov.(ff);
0930         case 'RETURN_RAW_DER'
0931             opt.opf.return_raw_der = ov.(ff);
0932         case 'FMC_ALG'
0933             opt.fmincon.alg = ov.(ff);
0934         case 'KNITRO_OPT'
0935             opt.knitro.opt = ov.(ff);
0936         case 'IPOPT_OPT'
0937             opt.ipopt.opt = ov.(ff);
0938         case 'MNS_FEASTOL'
0939             opt.minopf.feastol = ov.(ff);
0940         case 'MNS_ROWTOL'
0941             opt.minopf.rowtol = ov.(ff);
0942         case 'MNS_XTOL'
0943             opt.minopf.xtol = ov.(ff);
0944         case 'MNS_MAJDAMP'
0945             opt.minopf.majdamp = ov.(ff);
0946         case 'MNS_MINDAMP'
0947             opt.minopf.mindamp = ov.(ff);
0948         case 'MNS_PENALTY_PARM'
0949             opt.minopf.penalty = ov.(ff);
0950         case 'MNS_MAJOR_IT'
0951             opt.minopf.major_it = ov.(ff);
0952         case 'MNS_MINOR_IT'
0953             opt.minopf.minor_it = ov.(ff);
0954         case 'MNS_MAX_IT'
0955             opt.minopf.max_it = ov.(ff);
0956         case 'MNS_VERBOSITY'
0957             opt.minopf.verbosity = ov.(ff);
0958         case 'MNS_CORE'
0959             opt.minopf.core = ov.(ff);
0960         case 'MNS_SUPBASIC_LIM'
0961             opt.minopf.supbasic_lim = ov.(ff);
0962         case 'MNS_MULT_PRICE'
0963             opt.minopf.mult_price = ov.(ff);
0964         case 'FORCE_PC_EQ_P0'
0965             opt.sopf.force_Pc_eq_P0 = ov.(ff);
0966         case 'PDIPM_FEASTOL'
0967             opt.mips.feastol = ov.(ff);
0968             opt.pdipm.feastol = ov.(ff);
0969         case 'PDIPM_GRADTOL'
0970             opt.mips.gradtol = ov.(ff);
0971             opt.pdipm.gradtol = ov.(ff);
0972         case 'PDIPM_COMPTOL'
0973             opt.mips.comptol = ov.(ff);
0974             opt.pdipm.comptol = ov.(ff);
0975         case 'PDIPM_COSTTOL'
0976             opt.mips.costtol = ov.(ff);
0977             opt.pdipm.costtol = ov.(ff);
0978         case 'PDIPM_MAX_IT'
0979             opt.mips.max_it = ov.(ff);
0980             opt.pdipm.max_it = ov.(ff);
0981         case 'SCPDIPM_RED_IT'
0982             opt.mips.sc.red_it = ov.(ff);
0983             opt.pdipm.sc.red_it = ov.(ff);
0984         case 'TRALM_FEASTOL'
0985             opt.tralm.feastol = ov.(ff);
0986         case 'TRALM_PRIMETOL'
0987             opt.tralm.primaltol = ov.(ff);
0988         case 'TRALM_DUALTOL'
0989             opt.tralm.dualtol = ov.(ff);
0990         case 'TRALM_COSTTOL'
0991             opt.tralm.costtol = ov.(ff);
0992         case 'TRALM_MAJOR_IT'
0993             opt.tralm.major_it = ov.(ff);
0994         case 'TRALM_MINOR_IT'
0995             opt.tralm.minor_it = ov.(ff);
0996         case 'SMOOTHING_RATIO'
0997             opt.pdipm.sc.smooth_ratio = ov.(ff);
0998             opt.tralm.smooth_ratio = ov.(ff);
0999         case 'CPLEX_LPMETHOD'
1000             opt.cplex.lpmethod = ov.(ff);
1001         case 'CPLEX_QPMETHOD'
1002             opt.cplex.qpmethod = ov.(ff);
1003         case 'CPLEX_OPT'
1004             opt.cplex.opt = ov.(ff);
1005         case 'MOSEK_LP_ALG'
1006             opt.mosek.lp_alg = ov.(ff);
1007         case 'MOSEK_MAX_IT'
1008             opt.mosek.max_it = ov.(ff);
1009         case 'MOSEK_GAP_TOL'
1010             opt.mosek.gap_tol = ov.(ff);
1011         case 'MOSEK_MAX_TIME'
1012             opt.mosek.max_time = ov.(ff);
1013         case 'MOSEK_NUM_THREADS'
1014             opt.mosek.num_threads = ov.(ff);
1015         case 'MOSEK_OPT'
1016             opt.mosek.opt = ov.(ff);
1017         case 'GRB_METHOD'
1018             opt.gurobi.method = ov.(ff);
1019         case 'GRB_TIMELIMIT'
1020             opt.gurobi.timelimit = ov.(ff);
1021         case 'GRB_THREADS'
1022             opt.gurobi.threads = ov.(ff);
1023         case 'GRB_OPT'
1024             opt.gurobi.opt = ov.(ff);
1025         otherwise
1026             error('mpoption: ''%s'' is not a valid old-style option name', ff);
1027     end
1028 end
1029 % ov
1030 
1031 
1032 %%-------------------------------------------------------------------
1033 function opt_s = mpoption_v2s(opt_v)
1034 if DEBUG, fprintf('mpoption_v2s()\n'); end
1035 opt_s = mpoption_default();
1036 errstr = 'mpoption: %g is not a valid value for the old-style ''%s'' option';
1037 switch opt_v(1)                                 %% PF_ALG
1038     case 1
1039         opt_s.pf.alg = 'NR';        %% Newton's method
1040     case 2
1041         opt_s.pf.alg = 'FDXB';      %% fast-decoupled (XB version)
1042     case 3
1043         opt_s.pf.alg = 'FDBX';      %% fast-decoupled (BX version)
1044     case 4
1045         opt_s.pf.alg = 'GS';        %% Gauss-Seidel
1046     otherwise
1047         error(errstr, opt_v(1), 'PF_ALG');
1048 end
1049 opt_s.pf.tol                = opt_v(2);         %% PF_TOL
1050 opt_s.pf.nr.max_it          = opt_v(3);         %% PF_MAX_IT
1051 opt_s.pf.fd.max_it          = opt_v(4);         %% PF_MAX_IT_FD
1052 opt_s.pf.gs.max_it          = opt_v(5);         %% PF_MAX_IT_GS
1053 opt_s.pf.enforce_q_lims     = opt_v(6);         %% ENFORCE_Q_LIMS
1054 switch opt_v(10)                                %% PF_DC
1055     case 0
1056         opt_s.model = 'AC';
1057     case 1
1058         opt_s.model = 'DC';
1059     otherwise
1060         error(errstr, opt_v(10), 'PF_DC');
1061 end
1062 switch opt_v(11)                                %% OPF_ALG
1063     case 0
1064         opt_s.opf.ac.solver = 'DEFAULT';
1065     case 500
1066         opt_s.opf.ac.solver = 'MINOPF';
1067     case 520
1068         opt_s.opf.ac.solver = 'FMINCON';
1069     case {540, 545}
1070         opt_s.opf.ac.solver = 'PDIPM';
1071     case 550
1072         opt_s.opf.ac.solver = 'TRALM';
1073     case {560, 565}
1074         opt_s.opf.ac.solver = 'MIPS';
1075     case 580
1076         opt_s.opf.ac.solver = 'IPOPT';
1077     case 600
1078         opt_s.opf.ac.solver = 'KNITRO';
1079     otherwise
1080         error(errstr, opt_v(11), 'OPF_ALG');
1081 end
1082 opt_s.opf.violation         = opt_v(16);        %% OPF_VIOLATION
1083 
1084 opt_s.fmincon.tol_x         = opt_v(17);        %% CONSTR_TOL_X
1085 opt_s.fmincon.tol_f         = opt_v(18);        %% CONSTR_TOL_F
1086 opt_s.fmincon.max_it        = opt_v(19);        %% CONSTR_MAX_IT
1087 
1088 opt_s.knitro.tol_x          = opt_v(17);        %% CONSTR_TOL_X
1089 opt_s.knitro.tol_f          = opt_v(18);        %% CONSTR_TOL_F
1090 
1091 switch opt_v(24)                                %% OPF_FLOW_LIM
1092     case 0
1093         opt_s.opf.flow_lim = 'S';   %% apparent power (MVA)
1094     case 1
1095         opt_s.opf.flow_lim = 'P';   %% real power (MW)
1096     case 2
1097         opt_s.opf.flow_lim = 'I';   %% current magnitude (MVA @ 1 p.u. voltage)
1098     otherwise
1099         error(errstr, opt_v(10), 'PF_DC');
1100 end
1101 
1102 opt_s.opf.ignore_angle_lim  = opt_v(25);        %% OPF_IGNORE_ANG_LIM
1103 
1104 switch opt_v(26)                                %% OPF_ALG_DC
1105     case 0
1106         opt_s.opf.dc.solver = 'DEFAULT';
1107     case 100
1108         opt_s.opf.dc.solver = 'BPMPD';
1109     case {200, 250}
1110         opt_s.opf.dc.solver = 'MIPS';
1111     case 300
1112         opt_s.opf.dc.solver = 'OT';     %% QUADPROG, LINPROG
1113     case 400
1114         opt_s.opf.dc.solver = 'IPOPT';
1115     case 500
1116         opt_s.opf.dc.solver = 'CPLEX';
1117     case 600
1118         opt_s.opf.dc.solver = 'MOSEK';
1119     case 700
1120         opt_s.opf.dc.solver = 'GUROBI';
1121     otherwise
1122         error(errstr, opt_v(26), 'OPF_ALG_DC');
1123 end
1124 
1125 opt_s.verbose               = opt_v(31);        %% VERBOSE
1126 opt_s.out.all               = opt_v(32);        %% OUT_ALL
1127 opt_s.out.sys_sum           = opt_v(33);        %% OUT_SYS_SUM
1128 opt_s.out.area_sum          = opt_v(34);        %% OUT_AREA_SUM
1129 opt_s.out.bus               = opt_v(35);        %% OUT_BUS
1130 opt_s.out.branch            = opt_v(36);        %% OUT_BRANCH
1131 opt_s.out.gen               = opt_v(37);        %% OUT_GEN
1132 opt_s.out.lim.all           = opt_v(38);        %% OUT_ALL_LIM
1133 opt_s.out.lim.v             = opt_v(39);        %% OUT_V_LIM
1134 opt_s.out.lim.line          = opt_v(40);        %% OUT_LINE_LIM
1135 opt_s.out.lim.pg            = opt_v(41);        %% OUT_PG_LIM
1136 opt_s.out.lim.qg            = opt_v(42);        %% OUT_QG_LIM
1137 opt_s.out.force             = opt_v(44);        %% OUT_FORCE
1138 
1139 opt_s.opf.return_raw_der    = opt_v(52);        %% RETURN_RAW_DER
1140 
1141 opt_s.fmincon.alg           = opt_v(55);        %% FMC_ALG
1142 opt_s.knitro.opt            = opt_v(58);        %% KNITRO_OPT
1143 opt_s.ipopt.opt             = opt_v(60);        %% IPOPT_OPT
1144 
1145 opt_s.minopf.feastol        = opt_v(61);        %% MNS_FEASTOL
1146 opt_s.minopf.rowtol         = opt_v(62);        %% MNS_ROWTOL
1147 opt_s.minopf.xtol           = opt_v(63);        %% MNS_XTOL
1148 opt_s.minopf.majdamp        = opt_v(64);        %% MNS_MAJDAMP
1149 opt_s.minopf.mindamp        = opt_v(65);        %% MNS_MINDAMP
1150 opt_s.minopf.penalty        = opt_v(66);        %% MNS_PENALTY_PARM
1151 opt_s.minopf.major_it       = opt_v(67);        %% MNS_MAJOR_IT
1152 opt_s.minopf.minor_it       = opt_v(68);        %% MNS_MINOR_IT
1153 opt_s.minopf.max_it         = opt_v(69);        %% MNS_MAX_IT
1154 opt_s.minopf.verbosity      = opt_v(70);        %% MNS_VERBOSITY
1155 opt_s.minopf.core           = opt_v(71);        %% MNS_CORE
1156 opt_s.minopf.supbasic_lim   = opt_v(72);        %% MNS_SUPBASIC_LIM
1157 opt_s.minopf.mult_price     = opt_v(73);        %% MNS_MULT_PRICE
1158 
1159 opt_s.sopf.force_Pc_eq_P0   = opt_v(80);        %% FORCE_PC_EQ_P0, for c3sopf
1160 
1161 if (opt_v(11) == 565 && opt_v(10) == 0) || (opt_v(26) == 250 && opt_v(10) == 1)
1162     opt_s.mips.step_control = 1;
1163 end
1164 opt_s.mips.feastol          = opt_v(81);        %% PDIPM_FEASTOL
1165 opt_s.mips.gradtol          = opt_v(82);        %% PDIPM_GRADTOL
1166 opt_s.mips.comptol          = opt_v(83);        %% PDIPM_COMPTOL
1167 opt_s.mips.costtol          = opt_v(84);        %% PDIPM_COSTTOL
1168 opt_s.mips.max_it           = opt_v(85);        %% PDIPM_MAX_IT
1169 opt_s.mips.sc.red_it        = opt_v(86);        %% SCPDIPM_RED_IT
1170 
1171 opt_s.pdipm.feastol         = opt_v(81);        %% PDIPM_FEASTOL
1172 opt_s.pdipm.gradtol         = opt_v(82);        %% PDIPM_GRADTOL
1173 opt_s.pdipm.comptol         = opt_v(83);        %% PDIPM_COMPTOL
1174 opt_s.pdipm.costtol         = opt_v(84);        %% PDIPM_COSTTOL
1175 opt_s.pdipm.max_it          = opt_v(85);        %% PDIPM_MAX_IT
1176 opt_s.pdipm.sc.red_it       = opt_v(86);        %% SCPDIPM_RED_IT
1177 opt_s.pdipm.sc.smooth_ratio = opt_v(93);        %% SMOOTHING_RATIO
1178 if opt_v(11) == 545 && opt_v(10) == 0
1179     opt_s.pdipm.step_control = 1;
1180 end
1181 
1182 opt_s.tralm.feastol         = opt_v(87);        %% TRALM_FEASTOL
1183 opt_s.tralm.primaltol       = opt_v(88);        %% TRALM_PRIMETOL
1184 opt_s.tralm.dualtol         = opt_v(89);        %% TRALM_DUALTOL
1185 opt_s.tralm.costtol         = opt_v(90);        %% TRALM_COSTTOL
1186 opt_s.tralm.major_it        = opt_v(91);        %% TRALM_MAJOR_IT
1187 opt_s.tralm.minor_it        = opt_v(92);        %% TRALM_MINOR_IT
1188 opt_s.tralm.smooth_ratio    = opt_v(93);        %% SMOOTHING_RATIO
1189 
1190 opt_s.cplex.lpmethod        = opt_v(95);        %% CPLEX_LPMETHOD
1191 opt_s.cplex.qpmethod        = opt_v(96);        %% CPLEX_QPMETHOD
1192 opt_s.cplex.opt             = opt_v(97);        %% CPLEX_OPT
1193 
1194 opt_s.mosek.lp_alg          = opt_v(111);       %% MOSEK_LP_ALG
1195 opt_s.mosek.max_it          = opt_v(112);       %% MOSEK_MAX_IT
1196 opt_s.mosek.gap_tol         = opt_v(113);       %% MOSEK_GAP_TOL
1197 opt_s.mosek.max_time        = opt_v(114);       %% MOSEK_MAX_TIME
1198 opt_s.mosek.num_threads     = opt_v(115);       %% MOSEK_NUM_THREADS
1199 opt_s.mosek.opt             = opt_v(116);       %% MOSEK_OPT
1200 
1201 opt_s.gurobi.method         = opt_v(121);       %% GRB_METHOD
1202 opt_s.gurobi.timelimit      = opt_v(122);       %% GRB_TIMELIMIT
1203 opt_s.gurobi.threads        = opt_v(123);       %% GRB_THREADS
1204 opt_s.gurobi.opt            = opt_v(124);       %% GRB_OPT
1205 
1206 
1207 %%-------------------------------------------------------------------
1208 function opt_v = mpoption_s2v(opt_s)
1209 if DEBUG, fprintf('mpoption_s2v()\n'); end
1210 %% PF_ALG
1211 old = mpoption_old;
1212 switch upper(opt_s.pf.alg)
1213     case 'NR'
1214         PF_ALG = 1;
1215     case 'FDXB'
1216         PF_ALG = 2;
1217     case 'FDBX'
1218         PF_ALG = 3;
1219     case 'GS'
1220         PF_ALG = 4;
1221 end
1222 
1223 %% PF_DC
1224 if strcmp(upper(opt_s.model), 'DC')
1225     PF_DC = 1;
1226 else
1227     PF_DC = 0;
1228 end
1229 
1230 %% OPF_ALG
1231 switch upper(opt_s.opf.ac.solver)
1232     case 'DEFAULT'
1233         OPF_ALG = 0;
1234     case 'MINOPF'
1235         OPF_ALG = 500;
1236     case 'FMINCON'
1237         OPF_ALG = 520;
1238     case 'PDIPM'
1239         if isfield(opt_s, 'pdipm') && opt_s.pdipm.step_control
1240             OPF_ALG = 545;
1241         else
1242             OPF_ALG = 540;
1243         end
1244     case 'TRALM'
1245         OPF_ALG = 550;
1246     case 'MIPS'
1247         if opt_s.mips.step_control
1248             OPF_ALG = 565;
1249         else
1250             OPF_ALG = 560;
1251         end
1252     case 'IPOPT'
1253         OPF_ALG = 580;
1254     case 'KNITRO'
1255         OPF_ALG = 600;
1256 end
1257 
1258 %% FMINCON, Knitro tol_x, tol_f, max_it
1259 if strcmp(upper(opt_s.opf.ac.solver), 'KNITRO') && isfield(opt_s, 'knitro')
1260     CONSTR_TOL_X = opt_s.knitro.tol_x;
1261     CONSTR_TOL_F = opt_s.knitro.tol_f;
1262 elseif isfield(opt_s, 'fmincon')
1263     CONSTR_TOL_X  = opt_s.fmincon.tol_x;
1264     CONSTR_TOL_F  = opt_s.fmincon.tol_f;
1265 else
1266     CONSTR_TOL_X = old(17);
1267     CONSTR_TOL_F = old(18);
1268 end
1269 if isfield(opt_s, 'fmincon')
1270     CONSTR_MAX_IT   = opt_s.fmincon.max_it;
1271     FMC_ALG         = opt_s.fmincon.alg;
1272 else
1273     CONSTR_MAX_IT   = old(19);
1274     FMC_ALG         = old(55);
1275 end
1276 
1277 %% OPF_FLOW_LIM
1278 switch upper(opt_s.opf.flow_lim)
1279     case 'S'
1280         OPF_FLOW_LIM = 0;
1281     case {'P', '2'}
1282         OPF_FLOW_LIM = 1;
1283     case 'I'
1284         OPF_FLOW_LIM = 2;
1285 end
1286 
1287 %% OPF_ALG_DC
1288 switch upper(opt_s.opf.dc.solver)
1289     case 'DEFAULT'
1290         OPF_ALG_DC = 0;
1291     case 'BPMPD'
1292         OPF_ALG_DC = 100;
1293     case 'MIPS'
1294         if opt_s.mips.step_control
1295             OPF_ALG_DC = 250;
1296         else
1297             OPF_ALG_DC = 200;
1298         end
1299     case 'OT'
1300         OPF_ALG_DC = 300;
1301     case 'IPOPT'
1302         OPF_ALG_DC = 400;
1303     case 'CPLEX'
1304         OPF_ALG_DC = 500;
1305     case 'MOSEK'
1306         OPF_ALG_DC = 600;
1307     case 'GUROBI'
1308         OPF_ALG_DC = 700;
1309 end
1310 
1311 %% KNITRO_OPT
1312 if isfield(opt_s, 'knitro')
1313     KNITRO_OPT  = opt_s.knitro.opt;
1314 else
1315     KNITRO_OPT  = old(58);
1316 end
1317 
1318 %% IPOPT_OPT
1319 if isfield(opt_s, 'ipopt')
1320     IPOPT_OPT  = opt_s.ipopt.opt;
1321 else
1322     IPOPT_OPT  = old(58);
1323 end
1324 
1325 %% MINOPF options
1326 if isfield(opt_s, 'minopf')
1327     MINOPF_OPTS = [
1328         opt_s.minopf.feastol;   %% 61 - MNS_FEASTOL
1329         opt_s.minopf.rowtol;    %% 62 - MNS_ROWTOL
1330         opt_s.minopf.xtol;      %% 63 - MNS_XTOL
1331         opt_s.minopf.majdamp;   %% 64 - MNS_MAJDAMP
1332         opt_s.minopf.mindamp;   %% 65 - MNS_MINDAMP
1333         opt_s.minopf.penalty;   %% 66 - MNS_PENALTY_PARM
1334         opt_s.minopf.major_it;  %% 67 - MNS_MAJOR_IT
1335         opt_s.minopf.minor_it;  %% 68 - MNS_MINOR_IT
1336         opt_s.minopf.max_it;    %% 69 - MNS_MAX_IT
1337         opt_s.minopf.verbosity; %% 70 - MNS_VERBOSITY
1338         opt_s.minopf.core;      %% 71 - MNS_CORE
1339         opt_s.minopf.supbasic_lim;  %% 72 - MNS_SUPBASIC_LIM
1340         opt_s.minopf.mult_price;%% 73 - MNS_MULT_PRICE
1341     ];
1342 else
1343     MINOPF_OPTS = old(61:73);
1344 end
1345 
1346 %% FORCE_PC_EQ_P0
1347 if isfield(opt_s, 'sopf') && isfield(opt_s.sopf, 'force_Pc_eq_P0')
1348     FORCE_PC_EQ_P0 = opt_s.sopf.force_Pc_eq_P0;
1349 else
1350     FORCE_PC_EQ_P0 = 0;
1351 end
1352 
1353 %% PDIPM options
1354 if isfield(opt_s, 'pdipm')
1355     PDIPM_OPTS = [
1356         opt_s.pdipm.feastol;    %% 81 - PDIPM_FEASTOL
1357         opt_s.pdipm.gradtol;    %% 82 - PDIPM_GRADTOL
1358         opt_s.pdipm.comptol;    %% 83 - PDIPM_COMPTOL
1359         opt_s.pdipm.costtol;    %% 84 - PDIPM_COSTTOL
1360         opt_s.pdipm.max_it;     %% 85 - PDIPM_MAX_IT
1361         opt_s.pdipm.sc.red_it;  %% 86 - SCPDIPM_RED_IT
1362     ];
1363 else
1364     PDIPM_OPTS = old(81:86);
1365 end
1366 
1367 %% TRALM options
1368 if isfield(opt_s, 'tralm')
1369     TRALM_OPTS = [
1370         opt_s.tralm.feastol;    %% 87 - TRALM_FEASTOL
1371         opt_s.tralm.primaltol;  %% 88 - TRALM_PRIMETOL
1372         opt_s.tralm.dualtol;    %% 89 - TRALM_DUALTOL
1373         opt_s.tralm.costtol;    %% 90 - TRALM_COSTTOL
1374         opt_s.tralm.major_it;   %% 91 - TRALM_MAJOR_IT
1375         opt_s.tralm.minor_it;   %% 92 - TRALM_MINOR_IT
1376     ];
1377 else
1378     TRALM_OPTS = old(87:92);
1379 end
1380 
1381 %% SMOOTHING_RATIO
1382 if strcmp(upper(opt_s.opf.ac.solver), 'TRALM') && isfield(opt_s, 'tralm')
1383     SMOOTHING_RATIO = opt_s.tralm.smooth_ratio;
1384 elseif isfield(opt_s, 'pdipm')
1385     SMOOTHING_RATIO = opt_s.pdipm.sc.smooth_ratio;
1386 else
1387     SMOOTHING_RATIO = old(93);
1388 end
1389 
1390 %% CPLEX options
1391 if isfield(opt_s, 'cplex')
1392     CPLEX_OPTS = [
1393         opt_s.cplex.lpmethod;   %% 95 - CPLEX_LPMETHOD
1394         opt_s.cplex.qpmethod;   %% 96 - CPLEX_QPMETHOD
1395         opt_s.cplex.opt;        %% 97 - CPLEX_OPT
1396     ];
1397 else
1398     CPLEX_OPTS = old(95:97);
1399 end
1400 
1401 %% MOSEK options
1402 if isfield(opt_s, 'mosek')
1403     MOSEK_OPTS = [
1404         opt_s.mosek.lp_alg;     %% 111 - MOSEK_LP_ALG
1405         opt_s.mosek.max_it;     %% 112 - MOSEK_MAX_IT
1406         opt_s.mosek.gap_tol;    %% 113 - MOSEK_GAP_TOL
1407         opt_s.mosek.max_time;   %% 114 - MOSEK_MAX_TIME
1408         opt_s.mosek.num_threads;%% 115 - MOSEK_NUM_THREADS
1409         opt_s.mosek.opt;        %% 116 - MOSEK_OPT
1410     ];
1411 else
1412     MOSEK_OPTS = old(111:116);
1413 end
1414 
1415 %% Gurobi options
1416 if isfield(opt_s, 'gurobi')
1417     GUROBI_OPTS = [
1418         opt_s.gurobi.method;    %% 121 - GRB_METHOD
1419         opt_s.gurobi.timelimit; %% 122 - GRB_TIMELIMIT
1420         opt_s.gurobi.threads;   %% 123 - GRB_THREADS
1421         opt_s.gurobi.opt;       %% 124 - GRB_OPT
1422     ];
1423 else
1424     GUROBI_OPTS = old(121:124);
1425 end
1426 
1427 opt_v = [
1428         %% power flow options
1429         PF_ALG;                 %% 1  - PF_ALG
1430         opt_s.pf.tol;           %% 2  - PF_TOL
1431         opt_s.pf.nr.max_it;     %% 3  - PF_MAX_IT
1432         opt_s.pf.fd.max_it;     %% 4  - PF_MAX_IT_FD
1433         opt_s.pf.gs.max_it;     %% 5  - PF_MAX_IT_GS
1434         opt_s.pf.enforce_q_lims;%% 6  - ENFORCE_Q_LIMS
1435         0;                      %% 7  - RESERVED7
1436         0;                      %% 8  - RESERVED8
1437         0;                      %% 9  - RESERVED9
1438         PF_DC;                  %% 10 - PF_DC
1439         
1440         %% OPF options
1441         OPF_ALG;                %% 11 - OPF_ALG
1442         0;                      %% 12 - RESERVED12 (was OPF_ALG_POLY = 100)
1443         0;                      %% 13 - RESERVED13 (was OPF_ALG_PWL = 200)
1444         0;                      %% 14 - RESERVED14 (was OPF_POLY2PWL_PTS = 10)
1445         0;                      %% 15 - OPF_NEQ (removed)
1446         opt_s.opf.violation;    %% 16 - OPF_VIOLATION
1447         CONSTR_TOL_X;           %% 17 - CONSTR_TOL_X
1448         CONSTR_TOL_F;           %% 18 - CONSTR_TOL_F
1449         CONSTR_MAX_IT;          %% 19 - CONSTR_MAX_IT
1450         old(20);                %% 20 - LPC_TOL_GRAD (removed)
1451         old(21);                %% 21 - LPC_TOL_X (removed)
1452         old(22);                %% 22 - LPC_MAX_IT (removed)
1453         old(23);                %% 23 - LPC_MAX_RESTART (removed)
1454         OPF_FLOW_LIM;           %% 24 - OPF_FLOW_LIM
1455         opt_s.opf.ignore_angle_lim; %% 25 - OPF_IGNORE_ANG_LIM
1456         OPF_ALG_DC;             %% 26 - OPF_ALG_DC
1457         0;                      %% 27 - RESERVED27
1458         0;                      %% 28 - RESERVED28
1459         0;                      %% 29 - RESERVED29
1460         0;                      %% 30 - RESERVED30
1461         
1462         %% output options
1463         opt_s.verbose;          %% 31 - VERBOSE
1464         opt_s.out.all;          %% 32 - OUT_ALL
1465         opt_s.out.sys_sum;      %% 33 - OUT_SYS_SUM
1466         opt_s.out.area_sum;     %% 34 - OUT_AREA_SUM
1467         opt_s.out.bus;          %% 35 - OUT_BUS
1468         opt_s.out.branch;       %% 36 - OUT_BRANCH
1469         opt_s.out.gen;          %% 37 - OUT_GEN
1470         opt_s.out.lim.all;      %% 38 - OUT_ALL_LIM
1471         opt_s.out.lim.v;        %% 39 - OUT_V_LIM
1472         opt_s.out.lim.line;     %% 40 - OUT_LINE_LIM
1473         opt_s.out.lim.pg;       %% 41 - OUT_PG_LIM
1474         opt_s.out.lim.qg;       %% 42 - OUT_QG_LIM
1475         0;                      %% 43 - RESERVED43 (was OUT_RAW)
1476         opt_s.out.force;        %% 44 - OUT_FORCE
1477         0;                      %% 45 - RESERVED45
1478         0;                      %% 46 - RESERVED46
1479         0;                      %% 47 - RESERVED47
1480         0;                      %% 48 - RESERVED48
1481         0;                      %% 49 - RESERVED49
1482         0;                      %% 50 - RESERVED50
1483         
1484         %% other options
1485         old(51);                %% 51 - SPARSE_QP (removed)
1486         opt_s.opf.return_raw_der;   %% 52 - RETURN_RAW_DER
1487         0;                      %% 53 - RESERVED53
1488         0;                      %% 54 - RESERVED54
1489         FMC_ALG;                %% 55 - FMC_ALG
1490         0;                      %% 56 - RESERVED56
1491         0;                      %% 57 - RESERVED57
1492         KNITRO_OPT;             %% 58 - KNITRO_OPT
1493         0;                      %% 59 - RESERVED59
1494         IPOPT_OPT;              %% 60 - IPOPT_OPT
1495         
1496         %% MINOPF options
1497         MINOPF_OPTS;            %% 61-73 - MNS_FEASTOL-MNS_MULT_PRICE
1498         0;                      %% 74 - RESERVED74
1499         0;                      %% 75 - RESERVED75
1500         0;                      %% 76 - RESERVED76
1501         0;                      %% 77 - RESERVED77
1502         0;                      %% 78 - RESERVED78
1503         0;                      %% 79 - RESERVED79
1504         FORCE_PC_EQ_P0;         %% 80 - FORCE_PC_EQ_P0, for c3sopf
1505         
1506         %% MIPS, PDIPM, SC-PDIPM, and TRALM options
1507         PDIPM_OPTS;             %% 81-86 - PDIPM_FEASTOL-SCPDIPM_RED_IT
1508         TRALM_OPTS;             %% 87-92 - TRALM_FEASTOL-TRALM_MINOR_IT
1509         SMOOTHING_RATIO;        %% 93 - SMOOTHING_RATIO
1510         0;                      %% 94 - RESERVED94
1511         
1512         %% CPLEX options
1513         CPLEX_OPTS;             %% 95-97 - CPLEX_LPMETHOD-CPLEX_OPT
1514         0;                      %% 98 - RESERVED98
1515         0;                      %% 99 - RESERVED99
1516         0;                      %% 100 - RESERVED100
1517         0;                      %% 101 - RESERVED101
1518         0;                      %% 102 - RESERVED102
1519         0;                      %% 103 - RESERVED103
1520         0;                      %% 104 - RESERVED104
1521         0;                      %% 105 - RESERVED105
1522         0;                      %% 106 - RESERVED106
1523         0;                      %% 107 - RESERVED107
1524         0;                      %% 108 - RESERVED108
1525         0;                      %% 109 - RESERVED109
1526         0;                      %% 110 - RESERVED110
1527 
1528         %% MOSEK options
1529         MOSEK_OPTS;             %% 111-116 - MOSEK_LP_ALG-MOSEK_OPT
1530         0;                      %% 117 - RESERVED117
1531         0;                      %% 118 - RESERVED118
1532         0;                      %% 119 - RESERVED119
1533         0;                      %% 120 - RESERVED120
1534 
1535         %% Gurobi options
1536         GUROBI_OPTS;            %% 121-124 - GRB_METHOD-GRB_OPT
1537     ];
1538 
1539 
1540 %%-------------------------------------------------------------------
1541 function optt = mpoption_default()
1542 if DEBUG, fprintf('mpoption_default()\n'); end
1543 persistent opt;             %% cache this for speed
1544 if ~isstruct(opt)
1545     opt = struct(...
1546         'v',                    mpoption_version, ...   %% version
1547         'model',                'AC', ...
1548         'pf',                   struct(...
1549             'alg',                  'NR', ...
1550             'current_balance',      0, ...
1551             'v_cartesian',          0, ...
1552             'tol',                  1e-8, ...
1553             'nr',                   struct(...
1554                 'max_it',               10, ...
1555                 'lin_solver',            '' ), ...
1556             'fd',                   struct(...
1557                 'max_it',               30  ), ...
1558             'gs',                   struct(...
1559                 'max_it',               1000    ), ...
1560             'radial',               struct(...
1561                 'max_it',               20   , ...
1562                 'vcorr',                 0  ), ...
1563             'enforce_q_lims',       0   ), ...
1564         'cpf',                  struct(...
1565             'parameterization',     3, ...
1566             'stop_at',              'NOSE', ...     %% 'NOSE', <lam val>, 'FULL'
1567             'enforce_p_lims',       0, ...
1568             'enforce_q_lims',       0, ...
1569             'enforce_v_lims',       0, ...
1570             'enforce_flow_lims',    0, ...
1571             'step',                 0.05, ...
1572             'step_min',             1e-4, ...
1573             'step_max',             0.2, ...
1574             'adapt_step',           0, ...
1575             'adapt_step_damping',   0.7, ...
1576             'adapt_step_tol',       1e-3, ...
1577             'target_lam_tol',       1e-5, ...
1578             'nose_tol',             1e-5, ...
1579             'p_lims_tol',           0.01, ...
1580             'q_lims_tol',           0.01, ...
1581             'v_lims_tol',           1e-4, ...
1582             'flow_lims_tol',        0.01, ...
1583             'plot',                 struct(...
1584                 'level',                0, ...
1585                 'bus',                  []  ), ...
1586             'user_callback',        ''  ), ...
1587         'opf',                  struct(...
1588             'ac',                   struct(...
1589                 'solver',               'DEFAULT'   ), ...
1590             'dc',                   struct(...
1591                 'solver',               'DEFAULT'   ), ...
1592             'current_balance',      0, ...
1593             'v_cartesian',          0, ...
1594             'violation',            5e-6, ...
1595             'use_vg',               0, ...
1596             'flow_lim',             'S', ...
1597             'ignore_angle_lim',     0, ...
1598             'softlims',             struct(...
1599                 'default',               1  ), ...
1600             'init_from_mpc',        -1, ...
1601             'start',                0, ...
1602             'return_raw_der',       0   ), ...
1603         'verbose',              1, ...
1604         'out',                  struct(...
1605             'all',                  -1, ...
1606             'sys_sum',              1, ...
1607             'area_sum',             0, ...
1608             'bus',                  1, ...
1609             'branch',               1, ...
1610             'gen',                  0, ...
1611             'lim',                  struct(...
1612                 'all',                  -1, ...
1613                 'v',                    1, ...
1614                 'line',                 1, ...
1615                 'pg',                   1, ...
1616                 'qg',                   1   ), ...
1617             'force',                0, ...
1618             'suppress_detail',      -1  ), ...
1619         'mips',                 struct(...  %% see mpoption_info_mips() for optional fields
1620             'step_control',         0, ...
1621             'linsolver',            '', ...
1622             'feastol',              0, ...
1623             'gradtol',              1e-6, ...
1624             'comptol',              1e-6, ...
1625             'costtol',              1e-6, ...
1626             'max_it',               150, ...
1627             'sc',                   struct(...
1628                 'red_it',               20  )), ...
1629         'exp',                  struct(... %% experimental options
1630             'sys_wide_zip_loads',   struct(...
1631                 'pw',                   [], ...
1632                 'qw',                   []  )) ...
1633     );
1634     opt_pkgs = mpoption_optional_pkgs();
1635     for k = 1:length(opt_pkgs)
1636         fname = ['mpoption_info_' opt_pkgs{k}];
1637         if exist(fname, 'file') == 2
1638             opt = nested_struct_copy(opt, feval(fname, 'D'));
1639         end
1640     end
1641 end
1642 optt = opt;
1643 
1644 %%-------------------------------------------------------------------
1645 function optt = mpoption_optional_fields()
1646 if DEBUG, fprintf('mpoption_optional_fields()\n'); end
1647 persistent opt;         %% cache this for speed
1648 if ~isstruct(opt)
1649     opt_pkgs = mpoption_optional_pkgs();
1650     opt = struct;
1651     for k = 1:length(opt_pkgs)
1652         fname = ['mpoption_info_' opt_pkgs{k}];
1653         if exist(fname, 'file') == 2
1654             opt = nested_struct_copy(opt, feval(fname, 'V'));
1655         end
1656     end
1657 end
1658 optt = opt;
1659 
1660 %% globals
1661 %%-------------------------------------------------------------------
1662 function v = mpoption_version
1663 v = 21;     %% version number of MATPOWER options struct
1664             %% (must be incremented every time structure is updated)
1665             %% v1   - first version based on struct (MATPOWER 5.0b1)
1666             %% v2   - added 'linprog' and 'quadprog' fields
1667             %% v3   - (forgot to increment v) added 'out.suppress_detail'
1668             %%        field
1669             %% v4   - (forgot to increment v) MIPS 1.1, added optional
1670             %%        fields to 'mips' options: xi, sigma, z0, alpha_min,
1671             %%        rho_min, rho_max, mu_threshold and max_stepsize
1672             %% v5   - (forgot to increment v) added 'opf.init_from_mpc'
1673             %%        field (MATPOWER 5.0)
1674             %% v6   - added 'clp' field
1675             %% v7   - added 'intlinprog' field
1676             %% v8   - MIPS 1.2, added 'linsolver' field to
1677             %%        'mips' options
1678             %% v9   - added 'exp' for experimental fields, specifically
1679             %%        'sys_wide_zip_loads.pw', 'sys_wide_zip_loads.qw'
1680             %% v10  - added 'most' field
1681             %% v11  - added 'cpf' options 'adapt_step_damping',
1682             %%        'enforce_p_lims', 'enforce_q_lims', 'target_lam_tol'
1683             %%        'nose_tol', 'p_lims_tol' and 'q_lims_tol',
1684             %%        removed option 'cpf.user_callback_args'
1685             %% v12  - added option 'opf.use_vg'
1686             %% v13  - added 'pf.radial.max_it', 'pf.radial.vcorr'
1687             %% v14  - added 'pf.nr.lin_solver'
1688             %% v15  - added 'cpf.enforce_v_lims', 'cpf.enforce_flow_lims',
1689             %%        'cpf.v_lims_tol', and 'cpf.flow_lims_tol'
1690             %% v16  - added 'opf.start' (deprecated 'opf.init_from_mpc')
1691             %% v17  - added 'knitro.maxit'
1692             %% v18  - added 'opf.current_balance' and 'opf.v_cartesian'
1693             %% v19  - added 'opf.softlims.default'
1694             %% v20  - added 'pf.current_balance' and 'pf.v_cartesian'
1695             %% v21  - added 'osqp' field
1696 
1697 %%-------------------------------------------------------------------
1698 function db_level = DEBUG
1699 db_level = 0;
1700 
1701 %%-------------------------------------------------------------------
1702 function pkgs = mpoption_optional_pkgs()
1703 pkgs = {...
1704     'clp', 'cplex', 'fmincon', 'gurobi', 'glpk', 'intlinprog', 'ipopt', ...
1705     'knitro', 'linprog', 'minopf', 'most', 'mosek', 'osqp', 'quadprog', ...
1706     'sdp_pf', 'sopf', 'tspopf', 'yalmip' ...
1707 };

Generated on Fri 09-Oct-2020 11:21:31 by m2html © 2005