Implement OOP block structure description

CITATION.cff

@@ -10,5 +10,5 @@ authors:
     orcid: "https://orcid.org/0000-0002-1987-9268"
     affiliation: "McGill University"
 title: "decargroup/dkpy"
-version: v0.1.8
+version: v0.1.9
 url: "https://github.com/decargroup/dkpy"

README.rst

@@ -131,12 +131,12 @@ If you use this software in your research, please cite it as below or see
 
 .. code-block:: bibtex
 
-    @software{dahdah_dkpy_2024,
+    @software{dahdah_dkpy_2025,
         title={{decargroup/dkpy}},
         doi={10.5281/zenodo.14511244},
         url={https://github.com/decargroup/dkpy},
         publisher={Zenodo},
         author={Steven Dahdah and James Richard Forbes},
-        version = {{v0.1.8}},
-        year={2024},
+        version = {{v0.1.9}},
+        year={2025},
     }

doc/conf.py

@@ -9,8 +9,8 @@
 project = "dkpy"
 copyright = "2024, Steven Dahdah and James Richard Forbes"
 author = "Steven Dahdah and James Richard Forbes"
-version = "0.1.8"
-release = "0.1.8"
+version = "0.1.9"
+release = "0.1.9"
 
 # -- General configuration ---------------------------------------------------
 # https://www.sphinx-doc.org/en/master/usage/configuration.html#general-configuration

doc/dkpy.rst

@@ -69,6 +69,23 @@ selecting the order in :func:`DScaleFit.fit`.
 
    dkpy.DScaleFitSlicot
 
+Uncertainty Block Structure
+===========================
+
+Supported uncertainty block elements are provided below. The prefered method for 
+specifying the uncertainty block description is via :class:`UncertaintyBlock` objects.
+An additional method for specifying the uncertainty block description is the MATLAB 
+two column array format (see `MATLAB documentation <https://www.mathworks.com/help/robust/ref/mussv.html>`)
+for users that are more comfortable with this format.
+
+.. autosummary::
+   :toctree: _autosummary/
+
+   dkpy.UncertaintyBlock
+   dkpy.RealDiagonalBlock
+   dkpy.ComplexDiagonalBlock
+   dkpy.ComplexFullBlock
+
 Extending ``dkpy``
 ==================
 

examples/1_example_dk_iter_fixed_order.py

@@ -19,7 +19,12 @@ def example_dk_iter_fixed_order():
     )
 
     omega = np.logspace(-3, 3, 61)
-    block_structure = np.array([[1, 1], [1, 1], [2, 2]])
+    # block_structure = np.array([[1, 1], [1, 1], [2, 2]])  # Alternative MATLAB descr.
+    block_structure = [
+        dkpy.ComplexFullBlock(1, 1),
+        dkpy.ComplexFullBlock(1, 1),
+        dkpy.ComplexFullBlock(2, 2),
+    ]
     K, N, mu, iter_results, info = dk_iter.synthesize(
         eg["P"],
         eg["n_y"],

examples/2_example_dk_iter_list_order.py

@@ -23,7 +23,12 @@ def example_dk_iter_list_order():
     )
 
     omega = np.logspace(-3, 3, 61)
-    block_structure = np.array([[1, 1], [1, 1], [2, 2]])
+    # block_structure = np.array([[1, 1], [1, 1], [2, 2]])  # Alternative MATLAB descr.
+    block_structure = [
+        dkpy.ComplexFullBlock(1, 1),
+        dkpy.ComplexFullBlock(1, 1),
+        dkpy.ComplexFullBlock(2, 2),
+    ]
     K, N, mu, iter_results, info = dk_iter.synthesize(
         eg["P"],
         eg["n_y"],

examples/3_example_dk_iter_auto_order.py

@@ -44,7 +44,12 @@ def example_dk_iter_auto_order():
     )
 
     omega = np.logspace(-3, 3, 61)
-    block_structure = np.array([[1, 1], [1, 1], [2, 2]])
+    # block_structure = np.array([[1, 1], [1, 1], [2, 2]])  # Alternative MATLAB descr.
+    block_structure = [
+        dkpy.ComplexFullBlock(1, 1),
+        dkpy.ComplexFullBlock(1, 1),
+        dkpy.ComplexFullBlock(2, 2),
+    ]
     K, N, mu, iter_results, info = dk_iter.synthesize(
         eg["P"],
         eg["n_y"],

examples/4_example_dk_iter_interactive.py

@@ -37,7 +37,12 @@ def example_dk_iter_interactive():
     )
 
     omega = np.logspace(-3, 3, 61)
-    block_structure = np.array([[1, 1], [1, 1], [2, 2]])
+    # block_structure = np.array([[1, 1], [1, 1], [2, 2]])  # Alternative MATLAB descr.
+    block_structure = [
+        dkpy.ComplexFullBlock(1, 1),
+        dkpy.ComplexFullBlock(1, 1),
+        dkpy.ComplexFullBlock(2, 2),
+    ]
     K, N, mu, iter_results, info = dk_iter.synthesize(
         eg["P"],
         eg["n_y"],

examples/5_example_dk_iteration_custom.py

@@ -70,7 +70,12 @@ def example_dk_iter_custom():
     )
 
     omega = np.logspace(-3, 3, 61)
-    block_structure = np.array([[1, 1], [1, 1], [2, 2]])
+    # block_structure = np.array([[1, 1], [1, 1], [2, 2]])  # Alternative MATLAB descr.
+    block_structure = [
+        dkpy.ComplexFullBlock(1, 1),
+        dkpy.ComplexFullBlock(1, 1),
+        dkpy.ComplexFullBlock(2, 2),
+    ]
     K, N, mu, iter_results, info = dk_iter.synthesize(
         eg["P"],
         eg["n_y"],

pyproject.toml

@@ -4,7 +4,7 @@ build-backend = "hatchling.build"
 
 [project]
 name = "dkpy"
-version = "0.1.8"
+version = "0.1.9"
 dependencies = [
     "numpy>=1.21.0",
     "scipy>=1.7.0",

src/dkpy/__init__.py

@@ -4,4 +4,5 @@
 from .dk_iteration import *
 from .d_scale_fit import *
 from .structured_singular_value import *
+from .uncertainty_structure import *
 from .utilities import *

src/dkpy/d_scale_fit.py

@@ -6,7 +6,7 @@
 ]
 
 import abc
-from typing import Optional, Tuple, Union
+from typing import Optional, Tuple, Union, Sequence
 import warnings
 
 import control
@@ -15,6 +15,13 @@
 import slycot
 
 from . import utilities
+from .uncertainty_structure import (
+    UncertaintyBlock,
+    RealDiagonalBlock,
+    ComplexDiagonalBlock,
+    ComplexFullBlock,
+    _convert_matlab_block_structure,
+)
 
 
 class DScaleFit(metaclass=abc.ABCMeta):
@@ -26,7 +33,7 @@ def fit(
         omega: np.ndarray,
         D_omega: np.ndarray,
         order: Union[int, np.ndarray] = 0,
-        block_structure: Optional[np.ndarray] = None,
+        block_structure: Optional[Union[Sequence[UncertaintyBlock], np.ndarray]] = None,
     ) -> Tuple[control.StateSpace, control.StateSpace]:
         """Fit D-scale magnitudes.
 
@@ -39,10 +46,8 @@ def fit(
             dimension.
         order : Union[int, np.ndarray]
             Transfer function order to fit. Can be specified per-entry.
-        block_structure : np.ndarray
-            2D array with 2 columns and as many rows as uncertainty blocks
-            in Delta. The columns represent the number of rows and columns in
-            each uncertainty block. See [#mussv]_.
+        block_structure : Union[Sequence[UncertaintyBlock], np.ndarray]
+            Uncertainty block structure description.
 
         Returns
         -------
@@ -80,24 +85,24 @@ class DScaleFitSlicot(DScaleFit):
     ... )
     >>> D, D_inv = dkpy.DScaleFitSlicot().fit(omega, D_omega, 2, block_structure)
     >>> print(control.ss2tf(D[0, 0]))
-    <TransferFunction>: sys[3641]$indexed$converted
+    <TransferFunction>: sys[4066]$indexed$converted
     Inputs (1): ['u[0]']
     Outputs (1): ['y[0]']
     <BLANKLINE>
     <BLANKLINE>
-    0.0157 s^2 + 0.257 s + 0.1391
-    -----------------------------
-      s^2 + 0.9658 s + 0.01424
+    0.01571 s^2 + 0.257 s + 0.1391
+    ------------------------------
+       s^2 + 0.9657 s + 0.01424
     <BLANKLINE>
     >>> print(control.ss2tf(D[1, 1]))
-    <TransferFunction>: sys[3641]$indexed$converted
+    <TransferFunction>: sys[4066]$indexed$converted
     Inputs (1): ['u[1]']
     Outputs (1): ['y[1]']
     <BLANKLINE>
     <BLANKLINE>
     0.01573 s^2 + 0.2574 s + 0.1394
     -------------------------------
-       s^2 + 0.9663 s + 0.01425
+      s^2 + 0.9663 s + 0.01425
     <BLANKLINE>
     """
 
@@ -106,8 +111,10 @@ def fit(
         omega: np.ndarray,
         D_omega: np.ndarray,
         order: Union[int, np.ndarray] = 0,
-        block_structure: Optional[np.ndarray] = None,
+        block_structure: Optional[Union[Sequence[UncertaintyBlock], np.ndarray]] = None,
     ) -> Tuple[control.StateSpace, control.StateSpace]:
+        if isinstance(block_structure, np.ndarray):
+            block_structure = _convert_matlab_block_structure(block_structure)
         # Get mask
         if block_structure is None:
             mask = -1 * np.ones((D_omega.shape[0], D_omega.shape[1]), dtype=int)
@@ -155,18 +162,18 @@ def fit(
         return ss, ss_inv
 
 
-def _mask_from_block_structure(block_structure: np.ndarray) -> np.ndarray:
+def _mask_from_block_structure(
+    block_structure: Sequence[UncertaintyBlock],
+) -> np.ndarray:
     """Create a mask from a specified block structure.
 
     Entries known to be zero are set to 0. Entries known to be one are set to
     1. Entries to be fit numerically are set to -1.
 
     Parameters
     ----------
-    block_structure : np.ndarray
-        2D array with 2 columns and as many rows as uncertainty blocks
-        in Delta. The columns represent the number of rows and columns in
-        each uncertainty block. See [#mussv]_.
+    block_structure : Sequence[UncertaintyBlock]
+        Sequence of uncertainty block objects.
 
     Returns
     -------
@@ -178,19 +185,24 @@ def _mask_from_block_structure(block_structure: np.ndarray) -> np.ndarray:
     ----------
     .. [#mussv] https://www.mathworks.com/help/robust/ref/mussv.html
     """
+    num_blocks = len(block_structure)
     X_lst = []
-    for i in range(block_structure.shape[0]):
-        if block_structure[i, 0] <= 0:
+    for i in range(num_blocks):
+        # Uncertainty block
+        block = block_structure[i]
+        # Square uncertainty block condition
+        is_block_square = block.num_inputs == block.num_outputs
+        if isinstance(block, RealDiagonalBlock):
             raise NotImplementedError("Real perturbations are not yet supported.")
-        if block_structure[i, 1] <= 0:
+        if isinstance(block, ComplexDiagonalBlock):
             raise NotImplementedError("Diagonal perturbations are not yet supported.")
-        if block_structure[i, 0] != block_structure[i, 1]:
+        if isinstance(block, ComplexFullBlock) and (not is_block_square):
             raise NotImplementedError("Nonsquare perturbations are not yet supported.")
         # Set last scaling to identity
-        if i == block_structure.shape[0] - 1:
-            X_lst.append(np.eye(block_structure[i, 0], dtype=int))
+        if i == num_blocks - 1:
+            X_lst.append(np.eye(block.num_inputs, dtype=int))
         else:
-            X_lst.append(-1 * np.eye(block_structure[i, 0], dtype=int))
+            X_lst.append(-1 * np.eye(block.num_inputs, dtype=int))
     X = scipy.linalg.block_diag(*X_lst)
     return X