Wrapper Routines#
Since v1.0.10, following wrapper routines are introduced for easy access to other power flow and optimal power flow solvers. The wrapper routines are:
PYPOWER:
DCPF1,PFlow1,DCOPF1, andACOPF1gurobi_optimods:
OPF
Following packages are required to call them:
pypowergurobi_optimods
[1]:
import numpy as np
import ams
[2]:
sp = ams.load(ams.get_case('5bus/pjm5bus_demo.json'),
no_output=True)
PYPOWER#
Power Flow#
[3]:
sp.PFlow1.run()
Building system matrices
Parsing OModel for <PFlow1>
Evaluating OModel for <PFlow1>
Finalizing OModel for <PFlow1>
<PFlow1> converged in 0.0075 seconds, -1 iteration with PYPOWER.
PYPOWER Version 5.1.18, 10-Apr-2025 -- AC Power Flow (Newton)
Newton's method power flow converged in 3 iterations.
[3]:
True
You can switch to other power flow algorithms.
Note that it is recommended to use update to update the config values.
There is a known issue with this routine. Fast-Decoupled (XB version) pf_alg=2 and Fast-Decoupled (BX version) pf_alg=3 are not fully supported yet.
[4]:
sp.PFlow1.config.update(pf_alg=4)
[5]:
sp.PFlow1.run()
<PFlow1> converged in 0.0175 seconds, -1 iteration with PYPOWER.
PYPOWER Version 5.1.18, 10-Apr-2025 -- AC Power Flow (Gauss-Seidel)
Gauss-Seidel power flow converged in 69 iterations.
[5]:
True
After successful solving, the results are mapped back to routine variables.
[6]:
sp.PFlow1.pg.v
[6]:
array([1.02944468, 1. , 3.2349 , 4.6651 , 0.1 ])
[7]:
sp.PFlow1.qg.v
[7]:
array([ 1.67217961, -0.02705768, 0.95026549, -0.40953028, 1.31140954])
Below is the active power on line flow
[8]:
sp.PFlow1.plf.v
[8]:
array([ 1.34144848, 1.16636881, -2.84926578, -0.22728223, 0.00756141,
-1.80078873, 1.34144848])
A tcost object (an instance of ExpressionCalc) is created in the power flow routines for quick evaluation of the total generation cost.
[9]:
sp.PFlow1.tcost.v
[9]:
3.0854544680774185
Optimal Power Flow#
Similarly, both DCOPF1 and ACOPF1 are wrapped to solve optimal power flow.
[10]:
sp.DCOPF1.run()
Parsing OModel for <DCOPF1>
Evaluating OModel for <DCOPF1>
Finalizing OModel for <DCOPF1>
<DCOPF1> converged in 0.0137 seconds, 10 iterations with PYPOWER.
PYPOWER Version 5.1.18, 10-Apr-2025 -- DC Optimal Power Flow
Python Interior Point Solver - PIPS, Version 1.0, 07-Feb-2011
Converged!
[10]:
True
In PYPOWER, the c0 term (the constant coefficient in the generator cost function) is always included in the objective, regardless of the generator’s commitment status. This means tcosts and obj.v can be different when some generators are not committed.
See pypower/opf_costfcn.py for implementation details: https://rwl.github.io/PYPOWER/api/pypower.opf_costfcn-pysrc.html
[11]:
sp.DCOPF1.obj.e_str
[11]:
'sum(c2 * pg**2 + c1 * pg + c0)'
[12]:
sp.DCOPF1.obj.v
[12]:
0.9585953752140786
[13]:
sp.DCOPF1.tcost.e_str
[13]:
'sum(mul(c2, pg**2))+ sum(mul(c1, pg))+ sum(mul(ug, c0))'
[14]:
sp.DCOPF1.tcost.v
[14]:
0.9585953752140786
LMP at each bus is stored in variable pi
[15]:
sp.DCOPF1.pi.v
[15]:
array([0.00000776, 0.000001 , 0.00003 , 0.00001571, 0.00000917])
Kuhn-Tucker multiplier on upper and lower Pg limits are stored in mu1 and mu2 respectively.
[16]:
sp.DCOPF1.mu1.v
[16]:
array([0. , 0. , 0. , 0.0000342, 0. , 0. ,
0. ])
[17]:
sp.DCOPF1.mu2.v
[17]:
array([0., 0., 0., 0., 0., 0., 0.])
[18]:
sp.ACOPF1.run()
Parsing OModel for <ACOPF1>
Evaluating OModel for <ACOPF1>
Finalizing OModel for <ACOPF1>
PYPOWER Version 5.1.18, 10-Apr-2025 -- AC Optimal Power Flow
Python Interior Point Solver - PIPS, Version 1.0, 07-Feb-2011
<ACOPF1> converged in 0.1843 seconds, 15 iterations with PYPOWER.
Converged!
In ACOPF1, both active and reactive power prices are available.
[19]:
sp.ACOPF1.pi.v
[19]:
array([0.00000782, 0.000001 , 0.00003 , 0.00001561, 0.00000921])
[20]:
sp.ACOPF1.piq.v
[20]:
array([ 0.00000021, -0. , 0. , 0.0000001 , 0. ])
gurobi-optimods#
In this library, optimal power flow is modeled and solved using Gurobi.
https://gurobi-optimods.readthedocs.io/en/stable/mods/opf/opf.html
[21]:
help(sp.OPF.run)
Help on method run in module ams.routines.grbopt:
run(**kwargs) method of ams.routines.grbopt.OPF instance
Run the OPF routine using gurobi-optimods.
This method invokes `gurobi-optimods.opf.solve_opf` to solve the OPF problem.
Parameters
----------
- opftype : str
Type of OPF to solve (default: 'AC').
- branch_switching : bool
Enable branch switching (default: False).
- min_active_branches : float
Defines the minimum number of branches that must be turned on when
branch switching is active, i.e. the minimum number of turned on
branches is equal to ``numbranches * min_active_branches``. Has no
effect if ``branch_switching`` is set to False.
- use_mip_start : bool
Use MIP start (default: False).
- time_limit : float
Time limit for the solver (default: 0.0, no limit).
[22]:
sp.OPF.run(opftype='DC', branch_switching=True, time_limit=5)
Parsing OModel for <OPF>
Evaluating OModel for <OPF>
Finalizing OModel for <OPF>
Set parameter Username
Set parameter LicenseID to value 2617183
Set parameter TimeLimit to value 5
Set parameter OptimalityTol to value 0.0001
Academic license - for non-commercial use only - expires 2026-02-02
Building case data structures from dictionary.
Buses.
Bus 4 ID 3 is the reference bus.
sumloadPd 1000.0 numPload 3
sumloadQd 328.69
5 buses
Branches.
Numbranches: 7 active: 7
Generators.
Number of generators: 5
Number of buses with gens: 5
summaxPg 11430.0 summaxQg 10657.5
Generator cost vectors.
Running DCOPF formulation.
In bound_zs constraint, N=6
DCOPF model constructed (0.00s).
Statistics for model 'DC_Formulation_Model':
Problem type : MIP
Linear constraint matrix : 47 rows, 38 columns, 126 nonzeros
Variable types : 31 continuous, 7 integer (7 binary)
Matrix range : [1e-02, 2e+02]
Objective range : [1e+00, 1e+00]
Bounds range : [2e-01, 1e+02]
RHS range : [2e+00, 6e+00]
Using mip start with all branches kept on.
Gurobi Optimizer version 12.0.1 build v12.0.1rc0 (mac64[arm] - Darwin 23.6.0 23G93)
CPU model: Apple M3 Pro
Thread count: 12 physical cores, 12 logical processors, using up to 12 threads
Non-default parameters:
TimeLimit 5
OptimalityTol 0.0001
Optimize a model with 47 rows, 38 columns and 126 nonzeros
Model fingerprint: 0xbf6f625b
Variable types: 31 continuous, 7 integer (7 binary)
Coefficient statistics:
Matrix range [1e-02, 2e+02]
Objective range [1e+00, 1e+00]
Bounds range [2e-01, 1e+02]
RHS range [2e+00, 6e+00]
User MIP start produced solution with objective 0.958595 (0.01s)
Loaded user MIP start with objective 0.958595
Presolve removed 14 rows and 16 columns
Presolve time: 0.00s
Presolved: 33 rows, 22 columns, 111 nonzeros
Variable types: 15 continuous, 7 integer (7 binary)
Root relaxation: objective 6.802167e-01, 29 iterations, 0.00 seconds (0.00 work units)
Nodes | Current Node | Objective Bounds | Work
Expl Unexpl | Obj Depth IntInf | Incumbent BestBd Gap | It/Node Time
0 0 0.68022 0 2 0.95860 0.68022 29.0% - 0s
H 0 0 0.7341163 0.68022 7.34% - 0s
0 0 0.73412 0 2 0.73412 0.73412 0.00% - 0s
Cutting planes:
MIR: 3
Explored 1 nodes (33 simplex iterations) in 0.03 seconds (0.00 work units)
Thread count was 12 (of 12 available processors)
Solution count 2: 0.734116 0.958595
Optimal solution found (tolerance 1.00e-04)
Best objective 7.341162855941e-01, best bound 7.341162855940e-01, gap 0.0000%
Objective value = 0.7341162855940513.
Solution quality statistics for model 'DC_Formulation_Model' :
Maximum violation:
Bound : 0.00000000e+00
Constraint : 1.62758695e-13 (Pdef_3_0_4)
Integrality : 0.00000000e+00
<OPF> converged in 0.0588 seconds, -1 iteration with gurobi-optimods.
[22]:
True
[23]:
sp.OPF.pg.v
[23]:
array([0.92351428, 0.2 , 1. , 0.6 , 7.27648572])
As of v2.3.2, gurobi-optimods does not expose LMP
[24]:
sp.OPF.pi.v
[24]:
array([0., 0., 0., 0., 0.])
In this routine, a decision variable uld is introduced to store the branch_switching decision.
[25]:
sp.OPF.uld.v
[25]:
array([1., 1., 1., 1., 0., 1., 1.])
This variable values will also update the target model variables.
[26]:
sp.Line.u.v
[26]:
array([1., 1., 1., 1., 0., 1., 1.])