PowerModels-compatible OPF

The POWER_MODELS backend provides PowerModels-style network formulations through a Python-native modeling stack. It is useful for alternative AC/DC formulations, convex relaxations, reusable domain extensions, and large OPF problems that benefit from external optimization solvers.

Despite the name, this backend does not run Julia or call PowerModels.jl. Its formulation names and modeling concepts are aligned with PowerModels.jl, while pwrs builds and solves the models using Pyomo or CVXPY.

Install optional dependencies

The base pwrs installation is sufficient for the default MATPOWER/MIPS backend. The PowerModels-compatible backend loads its modeling frameworks lazily and requires packages appropriate to the selected formulation and solver.

When working from the source repository, install the complete optional stack with:

uv sync --group all

The main combinations are:

Workload	Modeling layer	Solver packages
Nonlinear AC formulations	Pyomo	CyIpopt and Ipopt libraries
LP/QP formulations	Pyomo	HiGHS (highspy), GLPK, or CyIpopt
SOCP formulations	CVXPY	Clarabel, SCS, or MOSEK
SDP formulations	CVXPY	SCS or MOSEK, depending on problem and installation

CyIpopt also requires an Ipopt installation and Pyomo’s PyNumero ASL interface. MOSEK requires a separately supplied license. See Setup with external optimizers for Ipopt and optional HSL setup.

Run the default ACP formulation

Select the backend explicitly. Its default formulation is ACP, and its default solver selection resolves a nonlinear ACP model to Ipopt.

import pwrs as mp

mpopt = mp.MatpowerConfig()
mpopt.opf.backend = "POWER_MODELS"
mpopt.opf.power_models.formulation = "ACP"
mpopt.opf.power_models.solver = "DEFAULT"
mpopt.verbose = 0
mpopt.out.all = 0

result = mp.runopf(mp.case30(), mpopt)
if not result.success:
    raise RuntimeError("PowerModels-compatible ACP did not converge")

print(result.f)
print(result.bus)
print(result.gen)

ACP is the exact nonlinear AC polar formulation and is the normal starting point when the result must be AC feasible. The backend maps its native model solution back to the same MATPOWER bus, generator, and branch result matrices as the default backend.

Choose a formulation

Set opf.power_models.formulation to one of the following names. A relaxation or approximation solves a different mathematical problem from ACP; a lower objective value does not mean it found a better AC-feasible dispatch.

Family	Formulations	Purpose
Exact nonlinear AC	ACP, ACR, ACT, IVR	AC OPF with polar, rectangular, lifted tangent, or current-voltage representations.
Active-power/DC	DCP, DCMP, NFA, DCPLL	DC, modified DC, network-flow, and DC-with-loss approximations.
AC approximations and relaxations	LPACC, BFA, SOCBF, SOCWR, QCRM, QCLS	Linearized, branch-flow, second-order-cone, and quadratic-convex models built with Pyomo.
Explicit conic relaxations	SOCWRCONIC, SOCBFCONIC, SDPWRM, SPARSESDPWRM	SOCP and SDP models built with CVXPY for conic solvers.

The active-power formulations DCP, DCMP, NFA, and DCPLL require the DC model. rundcopf() supplies that setting:

mpopt = mp.MatpowerConfig()
mpopt.opf.backend = "POWER_MODELS"
mpopt.opf.power_models.formulation = "DCP"
mpopt.opf.power_models.solver = "DEFAULT"
mpopt.verbose = 0
mpopt.out.all = 0

result = mp.rundcopf(mp.case30(), mpopt)

The other formulations require the AC model and should be called with runopf(). All PowerModels-compatible AC formulations currently enforce apparent-power branch ratings, so opf.flow_lim must remain "S".

When choosing among formulations:

• use ACP first for an AC-feasible operating point;

• use ACR or IVR for formulation studies or different numerical behavior;

• use DCP for a conventional active-power approximation;

• use NFA for a simpler network-flow relaxation;

• use the SOC, QC, or SDP formulations to compute relaxation bounds or study convex models, not as interchangeable replacements for an AC solution.

Automatic solver selection

With opf.power_models.solver = "DEFAULT", the backend classifies the built model after applying extensions and chooses an installed compatible solver:

Problem class	Default candidates
Linear Pyomo model	HiGHS, then GLPK, then Ipopt
Convex separable QP Pyomo model	HiGHS, then Ipopt
Nonlinear Pyomo model	Ipopt
CVXPY SOCP or SDP model	Clarabel, then SCS

An explicit solver request never silently falls back. For example:

mpopt.opf.power_models.solver = "HIGHS"
mpopt.opf.power_models.highs_options = {
    "primal_feasibility_tolerance": 1e-7,
}

Supported explicit names are HIGHS, GLPK, and IPOPT for Pyomo models, and CLARABEL, SCS, and MOSEK for CVXPY conic models. The backend rejects a solver that cannot handle the classified problem. For example, HiGHS cannot solve a nonlinear ACP model.

Solver-specific native options are stored in:

• opf.power_models.highs_options and glpk_options;

• ipopt.opts;

• opf.power_models.clarabel_options, scs_options, and mosek_options.

Configure Ipopt for large AC OPF

For medium and large ACP studies, an efficient sparse linear solver can have a much larger effect than small Python-level tuning changes:

import pwrs as mp
from pwrs.data.pglibopf import pglib_opf_case1354_pegase

mpopt = mp.MatpowerConfig()
mpopt.opf.backend = "POWER_MODELS"
mpopt.opf.power_models.formulation = "ACP"
mpopt.opf.power_models.solver = "IPOPT"
mpopt.ipopt.opts.update(
    {
        "linear_solver": "ma27",
        "hsllib": "libcoinhsl.so",  # use the installed platform-specific path
        "print_level": 0,
    }
)
mpopt.verbose = 0
mpopt.out.all = 0

result = mp.runopf(pglib_opf_case1354_pegase(), mpopt)
if not result.success:
    raise RuntimeError("Ipopt failed to solve ACP")

MA27 and other HSL solvers require a separately obtained Coin-HSL library. Do not copy the example library name blindly: use the path and filename from your installation, and confirm Ipopt reports the requested linear solver.

For performance comparisons, measure a warm and a cold end-to-end runopf() call. result.et includes case preparation, model construction, solver setup, the native solve, and result reconstruction. The finer timing breakdown in result.raw helps distinguish solver time from modeling overhead.

Inspect solver metadata

The backend records how a problem was classified and solved under result.raw["output"]:

if result.raw is not None:
    output = result.raw["output"]
    solver = output["solver"]

    print(output["formulation"])
    print(output["implementation"])
    print(solver["requested"], solver["selected"])
    print(solver["problem_class"])
    print(solver["fallbacks"])
    print(output["max_constraint_violation"])
    print(output["timings"])

Important fields include:

• formulation and implementation (PYOMO or CVXPY);

• requested and selected solver names;

• base and final problem class, including extension-driven reclassification;

• unavailable-solver fallbacks made during default selection;

• termination status and maximum constraint violation;

• network-build, model-build, solver, and result-mapping timings.

These fields are backend diagnostics rather than part of the stable structured result API. Code that must work with both OPF backends should depend on result.success, result.f, and the solved matrices instead.

Add reusable model extensions

Extensions run after the base formulation is built and before the model is classified for solver selection. An extension callback receives the public PowerModel protocol and returns None. The protocol exposes semantic quantities and named registration methods without requiring application code to depend on an internal context class.

Prefer the public domain helpers when they express the required operation. They work with both Pyomo and CVXPY formulations and preserve MATPOWER units and external component indexes.

import pwrs as mp
from pwrs.power_models import PowerModel, constraint_active_generation_sum


def add_reserve_floor(problem: PowerModel) -> None:
    constraint_active_generation_sum(
        problem,
        "reserve_floor",
        generator_rows=[0, 1],
        minimum_mw=80.0,
    )


mpopt = mp.MatpowerConfig()
mpopt.opf.backend = "POWER_MODELS"
mpopt.opf.power_models.formulation = "ACP"
mpopt.opf.power_models.extensions = (add_reserve_floor,)
mpopt.verbose = 0
mpopt.out.all = 0

result = mp.runopf(mp.case30(), mpopt)

Generator and branch selectors passed to the domain helpers are zero-based rows in the original MATPOWER matrices; bus selectors use BUS_I. Power quantities use MW/MVAr/MVA, angles use degrees, and voltage magnitudes use per unit. Available helpers cover generator active/reactive bounds and setpoints, voltage-magnitude bounds and setpoints, branch active/reactive/apparent-power limits, voltage-angle differences, weighted interface active flow, and aggregate active generation.

Callback names are used as extension result keys. Wrap a callback in PowerModelExtension(name, callback) when it needs a stable name different from its Python function name. Extension names and all variable, constraint, and objective-term group names must be unique identifiers.

Inspect capabilities before adding constraints

Formulations do not all expose the same physical quantities. For example, active-power-only formulations have no reactive generation, and NFA has no voltage-angle quantity. Reusable callbacks should query capabilities instead of assuming that a backend variable exists:

from pwrs.power_models import PowerModel


def add_reactive_floor(problem: PowerModel) -> None:
    if not problem.has_quantity("qg"):
        raise NotImplementedError(
            f"{problem.formulation} does not expose reactive generation"
        )

    problem.add_bounded_constraint(
        "first_generator_reactive_floor",
        problem.quantity("qg", 0),
        lower=-0.5,  # internal per-unit value
    )

The stable protocol surface consists of:

• network, formulation, and backend ("pyomo" or "cvxpy");

• available_quantities, has_quantity(), and quantity();

• the named bounded/equality constraint methods and add_constraint();

• add_variable() and add_objective_term();

• the read-only semantic variables, constraints, expressions, and objective_terms registries.

Domain helper arguments use public MATPOWER units. Values passed directly to quantity() and backend-native expressions use the formulation’s internal units, normally per unit for power.

Add backend-native variables, constraints, and objective terms

Use the backend-specific modeling package when a domain helper is not sufficient. The callback can branch on problem.backend while keeping its lifecycle and results common:

from pwrs.power_models import PowerModel


def add_shortfall_penalty(problem: PowerModel) -> None:
    if problem.backend == "pyomo":
        import pyomo.environ as pyo

        shortfall = problem.add_variable(
            "shortfall",
            pyo.Var(domain=pyo.NonNegativeReals),
        )
        balance = pyo.Constraint(
            expr=problem.quantity("pg", 0) + shortfall >= 2.0
        )
    else:
        import cvxpy as cp

        shortfall = problem.add_variable(
            "shortfall",
            cp.Variable(nonneg=True),
        )
        balance = problem.quantity("pg", 0) + shortfall >= 2.0

    problem.add_constraint(
        "shortfall_balance",
        balance,
        senses="lower",
    )
    problem.add_objective_term(
        "shortfall_penalty",
        10_000.0 * shortfall,
    )

add_constraint() accepts one native constraint or a tuple of constraints. Its optional senses= argument accepts "lower", "upper", "equality", "bounded", or "native", either once for a whole group or as one value per constraint. Pyomo directions are inferred when unambiguous. Specify directions for native CVXPY inequalities so that reported dual signs retain domain meaning. CVXPY additions must remain valid under its DCP rules.

Auxiliary variables do not change the MATPOWER-compatible result.x layout. Objective additions must be scalar and are included before final model classification and solver selection. Consequently, an extension can change a linear model into a QP or nonlinear problem and can change which default solver is selected.

Inspect structured extension results

Each named callback has a typed entry in result.extensions:

extension = result.extensions["add_shortfall_penalty"]

print(extension.build_time)
print(extension.variables["shortfall"])
print(extension.objective_terms["shortfall_penalty"])

constraint = extension.constraints["shortfall_balance"]
print(constraint.senses)
print(constraint.dual)
print(constraint.native_dual)
print(constraint.violation)
print(constraint.max_violation)

dual follows the registered semantic direction, while native_dual is the modeling package’s value. The same data remains in result.raw["extensions"] for mapping compatibility, but new code should use the structured result.extensions API.

For direct model manipulation, build and solve a formulation explicitly:

import pyomo.environ as pyo

from pwrs.power_models import build_acp_model, solve_acp_model

problem = build_acp_model(mp.case30())
problem.add_constraint(
    "minimum_generation",
    pyo.Constraint(expr=problem.quantity("pg", 0) >= 1.0),
    senses="lower",
)
result = solve_acp_model(problem, mpopt)
print(result.extensions["__built_model__"].constraints)

Direct additions made before solving are reported under the reserved "__built_model__" extension name. Register every native variable and constraint with the public methods: unregistered active constraints are rejected to prevent a solve whose custom behavior is missing from diagnostics. See the API Reference for the domain helpers, formulation builders, and matching solve functions.

Current input and option limitations

The PowerModels-compatible backend deliberately rejects unsupported input instead of silently dropping it. At present:

• input must be a MATPOWER case; native PowerModels component dictionaries, multinetwork data, and multiconductor data are not supported;

• zero-impedance branches and negative branch thermal ratings are rejected;

• every active generator needs an active-power gencost row, using a polynomial or piecewise-linear model;

• nonempty MATPOWER generalized constraint/cost fields (A, l, u, N, H, Cw, fparm, z0, zl, or zu) are rejected;

• MATPOWER userfcn callbacks are rejected; use opf.power_models.extensions instead;

• only options compatible with the selected formulation are accepted, and AC formulations currently support apparent-power branch limits only;

• some nonstandard component types are not yet supported.

DC lines are supported subject to their bounds, loss, and cost validation. If a case depends on a feature in this list, use the MATPOWER-compatible backend or reformulate the study explicitly.

Troubleshooting

When a PowerModels-compatible solve fails:

1. Confirm the modeling package and requested solver are installed in the same environment used by uv run.

2. Check that the formulation uses the correct AC or DC model.

3. Leave the solver at DEFAULT to test capability-based selection.

4. Inspect solver fallback and problem-class metadata when a result exists.

5. For Ipopt, verify PyNumero ASL, CyIpopt, Ipopt, and any requested HSL library independently.

6. For a relaxation, distinguish solver infeasibility from the fact that its solution may not map to an AC-feasible operating point.

7. Remove user extensions and solve the base formulation before debugging an added constraint.

8. Compare the same case with ACP and the MATPOWER backend to separate an input problem from a formulation- or solver-specific issue.