Skip to content

Commit

Permalink
Merge pull request #80 from CURENT/misc
Browse files Browse the repository at this point in the history
Add static analysis matrices calculation
  • Loading branch information
jinningwang authored Apr 21, 2024
2 parents 2ff1ae9 + 4d7ca5c commit cab0acd
Show file tree
Hide file tree
Showing 12 changed files with 599 additions and 199 deletions.
427 changes: 332 additions & 95 deletions ams/core/matprocessor.py

Large diffs are not rendered by default.

2 changes: 1 addition & 1 deletion ams/interop/andes.py
Original file line number Diff line number Diff line change
Expand Up @@ -769,7 +769,7 @@ def receive(self, adsys=None, routine=None, no_update=False):
if no_update and (len(pname_to_update) > 0):
logger.info(f'Please update <{rtn.class_name}> parameters: {pname_to_update}')
elif len(pname_to_update) > 0:
rtn.update(params=pname_to_update, mat_make=False)
rtn.update(params=pname_to_update, build_mats=False)
return True


Expand Down
2 changes: 1 addition & 1 deletion ams/pypower/routines/pflow.py
Original file line number Diff line number Diff line change
Expand Up @@ -342,7 +342,7 @@ def dcpf(B, Pbus, Va0, ref, pv, pq):
Ray Zimmerman (PSERC Cornell)
"""
pvpq = np.matrix(np.r_[pv, pq])
pvpq = np.r_[pv, pq][np.newaxis, :]

# initialize result vector
Va = np.copy(Va0)
Expand Down
3 changes: 3 additions & 0 deletions ams/routines/dcopf.py
Original file line number Diff line number Diff line change
Expand Up @@ -315,6 +315,7 @@ def dc2ac(self, kloss=1.0, **kwargs):
kloss : float, optional
The loss factor for the conversion. Defaults to 1.2.
"""
exec_time = self.exec_time
if self.exec_time == 0 or self.exit_code != 0:
logger.warning(f'{self.class_name} is not executed successfully, quit conversion.')
return False

Check warning on line 321 in ams/routines/dcopf.py

View check run for this annotation

Codecov / codecov/patch

ams/routines/dcopf.py#L318-L321

Added lines #L318 - L321 were not covered by tests
Expand All @@ -330,6 +331,7 @@ def dc2ac(self, kloss=1.0, **kwargs):
ACOPF = self.system.ACOPF

Check warning on line 331 in ams/routines/dcopf.py

View check run for this annotation

Codecov / codecov/patch

ams/routines/dcopf.py#L331

Added line #L331 was not covered by tests
# run ACOPF
ACOPF.run()

Check warning on line 333 in ams/routines/dcopf.py

View check run for this annotation

Codecov / codecov/patch

ams/routines/dcopf.py#L333

Added line #L333 was not covered by tests
# self.exec_time += ACOPF.exec_time
# scale load back
self.system.StaticLoad.set(src='p0', attr='v', idx=pq_idx, value=pd0)
self.system.StaticLoad.set(src='q0', attr='v', idx=pq_idx, value=qd0)
Expand All @@ -349,6 +351,7 @@ def dc2ac(self, kloss=1.0, **kwargs):
self.vBus.parse()
self.vBus.optz.value = ACOPF.vBus.v
self.aBus.optz.value = ACOPF.aBus.v
self.exec_time = exec_time

Check warning on line 354 in ams/routines/dcopf.py

View check run for this annotation

Codecov / codecov/patch

ams/routines/dcopf.py#L351-L354

Added lines #L351 - L354 were not covered by tests

# --- set status ---
self.system.recent = self
Expand Down
9 changes: 9 additions & 0 deletions ams/routines/dcpf.py
Original file line number Diff line number Diff line change
Expand Up @@ -64,6 +64,12 @@ def __init__(self, system, config):
name='aBus', tex_name=r'a_{Bus}',
model='Bus', src='a',)

# --- line flow ---
self.plf = Var(info='Line flow',
unit='p.u.',
name='plf', tex_name=r'p_{lf}',
model='Line',)

def unpack(self, res):
"""
Unpack results from PYPOWER.
Expand All @@ -84,6 +90,9 @@ def unpack(self, res):
system.Slack.p.v = res['gen'][:system.Slack.n, 1] / mva # active power
system.Slack.q.v = res['gen'][:system.Slack.n, 2] / mva # reactive power

# --- Line ---
self.plf.optz.value = res['branch'][:, 13] / mva # line flow

# --- copy results from system algeb into routine algeb ---
for vname, var in self.vars.items():
owner = getattr(system, var.model) # instance of owner, Model or Group
Expand Down
20 changes: 10 additions & 10 deletions ams/routines/routine.py
Original file line number Diff line number Diff line change
Expand Up @@ -107,7 +107,7 @@ def get(self, src: str, idx, attr: str = 'v',
if isinstance(idx, (str, int)):
idx_u = [idx]
is_format = True
elif isinstance(idx, np.ndarray):
elif isinstance(idx, (np.ndarray, pd.Series)):
idx_u = idx.tolist()
elif isinstance(idx, list):
idx_u = idx.copy()
Expand Down Expand Up @@ -267,13 +267,13 @@ def init(self, force=False, no_code=True, **kwargs):

# --- force initialization ---
if force:
self.system.mats.make()
self.system.mats.build()

Check warning on line 270 in ams/routines/routine.py

View check run for this annotation

Codecov / codecov/patch

ams/routines/routine.py#L270

Added line #L270 was not covered by tests
for constr in self.constrs.values():
constr.is_disabled = False

# --- matrix build ---
if not self.system.mats.initialized:
self.system.mats.make()
self.system.mats.build()

# --- constraint check ---
_ = self._get_off_constrs()
Expand Down Expand Up @@ -492,24 +492,24 @@ def _register_attribute(self, key, value):
elif isinstance(value, RBaseService):
self.services[key] = value

def update(self, params=None, mat_make=True,):
def update(self, params=None, build_mats=True,):
"""
Update the values of Parameters in the optimization model.
This method is particularly important when some `RParams` are
linked with system matrices.
In such cases, setting `mat_make=True` is necessary to rebuild
In such cases, setting `build_mats=True` is necessary to rebuild
these matrices for the changes to take effect.
This is common in scenarios involving topology changes, connection statuses,
or load value modifications.
If unsure, it is advisable to use `mat_make=True` as a precautionary measure.
If unsure, it is advisable to use `build_mats=True` as a precautionary measure.
Parameters
----------
params: Parameter, str, or list
Parameter, Parameter name, or a list of parameter names to be updated.
If None, all parameters will be updated.
mat_make: bool
build_mats: bool
True to rebuild the system matrices. Set to False to speed up the process
if no system matrices are changed.
"""
Expand All @@ -519,7 +519,7 @@ def update(self, params=None, mat_make=True,):
sparams = []
if params is None:
sparams = [val for val in self.params.values()]
mat_make = True
build_mats = True
elif isinstance(params, Param):
sparams = [params]
elif isinstance(params, str):
Expand All @@ -532,8 +532,8 @@ def update(self, params=None, mat_make=True,):
if param.optz is None: # means no_parse=True
re_init = True
break
if mat_make:
self.system.mats.make()
if build_mats:
self.system.mats.build()
if re_init:
logger.warning(f"<{self.class_name}> reinit OModel due to non-parametric change.")
self.om.parsed = False
Expand Down
2 changes: 2 additions & 0 deletions ams/routines/rted.py
Original file line number Diff line number Diff line change
Expand Up @@ -184,6 +184,7 @@ def dc2ac(self, kloss=1.0, **kwargs):
kloss : float, optional
The loss factor for the conversion. Defaults to 1.2.
"""
exec_time = self.exec_time
if self.exec_time == 0 or self.exit_code != 0:
logger.warning(f'{self.class_name} is not executed successfully, quit conversion.')

Check warning on line 189 in ams/routines/rted.py

View check run for this annotation

Codecov / codecov/patch

ams/routines/rted.py#L189

Added line #L189 was not covered by tests
return False
Expand Down Expand Up @@ -228,6 +229,7 @@ def dc2ac(self, kloss=1.0, **kwargs):
self.vBus.parse()
self.vBus.optz.value = ACOPF.vBus.v
self.aBus.optz.value = ACOPF.aBus.v
self.exec_time = exec_time

# reset pmin, pmax, p0
self.system.StaticGen.set(src='pmin', attr='v', idx=pr_idx, value=pmin0)
Expand Down
9 changes: 9 additions & 0 deletions docs/source/release-notes.rst
Original file line number Diff line number Diff line change
Expand Up @@ -12,13 +12,22 @@ Pre-v1.0.0
v0.9.6 (2024-xx-xx)
-------------------

This patch release refactor and improve `MatProcessor`, where it support PTDF, LODF,
and OTDF for static analysis.

The reference can be found online "PowerWorld > Web Help > Sensitivities > Line
Outage Distribution Factors".

- Refactor DCPF, PFlow, and ACOPF
- Add a loss factor in ``RTED.dc2ac()``
- Add ``DCOPF.dc2ac()``
- Fix OModel parse status to ensure no_parsed params can be updated
- Fix and rerun ex2
- Format ``Routine.get()`` return type to be consistent with input idx type
- Remove unused ``Routine.prepare()``
- Refactor `MatProcessor` to separate matrix building
- Add Var `plf` in `DCPF`, `PFlow`, and `ACOPF` to store the line flow
- Add `build_ptdf`, `build_lodf`, and `build_otdf`

v0.9.5 (2024-03-25)
-------------------
Expand Down
96 changes: 47 additions & 49 deletions examples/demonstration/demo_AGC.ipynb

Large diffs are not rendered by default.

183 changes: 183 additions & 0 deletions tests/test_mats.py
Original file line number Diff line number Diff line change
@@ -0,0 +1,183 @@
import unittest
import numpy as np

from scipy.sparse import csr_matrix as c_sparse
from scipy.sparse import lil_matrix as l_sparse

import ams
from ams.core.matprocessor import MatProcessor, MParam


class TestMatProcessor(unittest.TestCase):
"""
Test functionality of MatProcessor.
"""

def setUp(self) -> None:
# cases is for testing PTDF, LODF, etc.
self.cases = ['matpower/case14.m',
'matpower/case39.m',
'matpower/case118.m']
self.ss = ams.load(ams.get_case("matpower/case300.m"),
default_config=True, no_output=True)
self.nR = self.ss.Region.n
self.nb = self.ss.Bus.n
self.nl = self.ss.Line.n
self.ng = self.ss.StaticGen.n
self.nsh = self.ss.Shunt.n
self.nD = self.ss.StaticLoad.n

self.mats = MatProcessor(self.ss)
self.mats.build()

def test_MParam(self):
"""
Test `MParam`.
"""
one_vec = MParam(v=c_sparse(np.ones(self.ss.Bus.n)))
# check if `_v` is `c_sparse` instance
self.assertIsInstance(one_vec._v, c_sparse)
# check if `v` is 1D-array
self.assertEqual(one_vec.v.shape, (self.ss.Bus.n,))

def test_cg(self):
"""
Test `Cg`.
"""
self.assertIsInstance(self.mats.Cg._v, (c_sparse, l_sparse))
self.assertIsInstance(self.mats.Cg.v, np.ndarray)
self.assertEqual(self.mats.Cg._v.max(), 1)
np.testing.assert_equal(self.mats.Cg._v.sum(axis=0), np.ones((1, self.ng)))

def test_cl(self):
"""
Test `Cl`.
"""
self.assertIsInstance(self.mats.Cl._v, (c_sparse, l_sparse))
self.assertIsInstance(self.mats.Cl.v, np.ndarray)
self.assertEqual(self.mats.Cl._v.max(), 1)
np.testing.assert_equal(self.mats.Cl._v.sum(axis=0), np.ones((1, self.nD)))

def test_csh(self):
"""
Test `Csh`.
"""
self.assertIsInstance(self.mats.Csh._v, (c_sparse, l_sparse))
self.assertIsInstance(self.mats.Csh.v, np.ndarray)
self.assertEqual(self.mats.Csh._v.max(), 1)
np.testing.assert_equal(self.mats.Csh._v.sum(axis=0), np.ones((1, self.nsh)))

def test_cft(self):
"""
Test `Cft`.
"""
self.assertIsInstance(self.mats.Cft._v, (c_sparse, l_sparse))
self.assertIsInstance(self.mats.Cft.v, np.ndarray)
self.assertEqual(self.mats.Cft._v.max(), 1)
np.testing.assert_equal(self.mats.Cft._v.sum(axis=0), np.zeros((1, self.nl)))

def test_b_calc(self):
"""
Test `b_calc`.
"""
b = self.mats._calc_b()
self.assertIsInstance(b, np.ndarray)
self.assertEqual(b.shape, (self.nl,))

def test_bf(self):
"""
Test `Bf`.
"""
self.assertIsInstance(self.mats.Bf._v, (c_sparse, l_sparse))
np.testing.assert_equal(self.mats.Bf._v.shape, (self.nl, self.nb))

def test_bbus(self):
"""
Test `Bbus`.
"""
self.assertIsInstance(self.mats.Bbus._v, (c_sparse, l_sparse))
np.testing.assert_equal(self.mats.Bbus._v.shape, (self.nb, self.nb))

def test_pfinj(self):
"""
Test `Pfinj`.
"""
self.assertIsInstance(self.mats.Pfinj._v, np.ndarray)
np.testing.assert_equal(self.mats.Pfinj._v.shape, (self.nl,))

def test_pbusinj(self):
"""
Test `Pbusinj`.
"""
self.assertIsInstance(self.mats.Pbusinj._v, np.ndarray)
np.testing.assert_equal(self.mats.Pbusinj._v.shape, (self.nb,))

def test_ptdf(self):
"""
Test `PTDF`.
"""

for case in self.cases:
ss = ams.load(ams.get_case(case),
setup=True, default_config=True, no_output=True)
ss.DCOPF.run(solver='ECOS')

ptdf = ss.mats.build_ptdf()

plf = ss.DCOPF.plf.v
plfc = ptdf@(ss.mats.Cg._v@ss.DCOPF.pg.v - ss.mats.Cl._v@ss.DCOPF.pd.v)
np.testing.assert_allclose(plf, plfc, atol=1e-2)

def test_lodf(self):
"""
Test `LODF`.
"""
for case in self.cases:
ss = ams.load(ams.get_case(case),
setup=True, default_config=True, no_output=True)
# build matrices
ss.mats.build()
_ = ss.mats.build_ptdf()
lodf = ss.mats.build_lodf()

# outage line
oline_idx = ss.Line.idx.v[1]
oline = ss.Line.idx2uid(oline_idx)

# pre-outage
ss.DCPF.run()
plf0 = ss.DCPF.plf.v.copy()

# post-outage
ss.Line.set(src='u', attr='v', idx=oline_idx, value=0)
ss.DCPF.update()
ss.DCPF.run()
plf1 = ss.DCPF.plf.v.copy()

dplf = plf1 - plf0
dplfc = -lodf[:, oline] * dplf[oline]

np.testing.assert_allclose(dplf, dplfc, atol=1e-7)

def test_otdf(self):
"""
Test `OTDF`.
"""

for case in self.cases:
ss = ams.load(ams.get_case(case),
setup=True, default_config=True, no_output=True)
# build matrices
ss.mats.build()

oline_idx = ss.Line.idx.v[1]

otdf = ss.mats.build_otdf(line=oline_idx)

ss.Line.set(src='u', attr='v', idx=oline_idx, value=0)
ss.DCPF.run()

plf = ss.DCPF.plf.v
plfc = otdf@(ss.mats.Cg._v@ss.DCPF.pg.v - ss.mats.Cl._v@ss.DCPF.pd.v)

np.testing.assert_allclose(plf, plfc, atol=1e-7)
1 change: 1 addition & 0 deletions tests/test_routine.py
Original file line number Diff line number Diff line change
Expand Up @@ -158,6 +158,7 @@ def test_dc2ac(self):
"""
self.ss.RTED.dc2ac()
self.assertTrue(self.ss.RTED.converted, "AC conversion failed!")
self.assertTrue(self.ss.RTED.exec_time > 0, "Execution time is not greater than 0.")

stg_idx = self.ss.StaticGen.get_idx()
pg_rted = self.ss.RTED.get(src='pg', attr='v', idx=stg_idx)
Expand Down
Loading

0 comments on commit cab0acd

Please sign in to comment.