ellalgo.oracles namespace

Submodules

ellalgo.oracles.ldlt_mgr module

LDLTMgr (LDLT Manager)

This code defines a class called LDLTMgr, which stands for LDLT Manager. The purpose of this class is to perform a special kind of matrix factorization called LDLT factorization on symmetric matrices. This factorization is useful in various mathematical and engineering applications, especially when dealing with linear systems and eigenvalue problems.

The main input for this class is a symmetric matrix, which can be provided either as a NumPy array or through a function that returns individual matrix elements. The primary output is a boolean value indicating whether the input matrix is positive definite (a special property of matrices). Additionally, the class can produce other outputs like a witness vector (if the matrix is not positive definite) and a square root of the matrix (if it is positive definite).

The LDLTMgr class achieves its purpose through several methods. The main method is ‘factor’, which performs the actual LDLT factorization. It does this by going through the matrix elements row by row, calculating and storing values in a special way. This process helps determine if the matrix is positive definite and allows for efficient calculations later.

An important aspect of this code is its use of “lazy evaluation”. This means it doesn’t need the entire matrix upfront but can work with just a function that provides matrix elements as needed. This can be more efficient for large matrices or when the matrix is defined by a formula rather than stored values.

The class also includes methods to check if the matrix is positive definite (‘is_spd’), calculate a witness vector if it’s not (‘witness’), and compute a symmetric quadratic form (‘sym_quad’). These additional functionalities make the class versatile for various matrix-related tasks.

One of the key transformations happening in this code is the factorization itself. It’s breaking down the original matrix into simpler parts (L, D, and L transpose), which can be used to solve problems more easily. This factorization is done in a way that avoids using square roots, which can be computationally expensive.

Overall, the LDLTMgr class provides a set of tools for working with symmetric matrices, with a focus on determining their properties (like positive definiteness) and performing useful calculations efficiently. It’s designed to be flexible and can handle both standard matrices and those defined by functions, making it useful in a variety of mathematical and engineering contexts.

class ellalgo.oracles.ldlt_mgr.LDLTMgr(ndim: int)

Bases: object

LDLT factorization (mainly for LMI oracles)

LDLTMgr is a class that performs the LDLT factorization for a given symmetric matrix. The LDLT factorization decomposes a symmetric matrix A into the product of a lower triangular matrix L, a diagonal matrix D, and the transpose of L. This factorization is useful for solving linear systems and eigenvalue problems. The class provides methods to perform the factorization, check if the matrix is positive definite, calculate a witness vector if it is not positive definite, and calculate the symmetric quadratic form.

• LDL^T square-root-free version

• Option allow semidefinite

• Cholesky-Banachiewicz style, row-based

• Lazy evaluation

• A matrix A in R^{m x m} is positive definite

  iff v^T A v > 0 for all v in R^n.

• O(p^3) per iteration, independent of ndim

Examples

>>> A = np.array([[1.0, 0.5, 0.5], [0.5, 1.25, 0.75], [0.5, 0.75, 1.5]])
>>> ldl = LDLTMgr(3)
>>> ldl.factorize(A)
True

factor(get_elem: Callable[[int, int], float]) → bool

Performs LDLT factorization using lazy element access.

The factorization proceeds row by row, computing diagonal entries and off-diagonal multipliers. If any diagonal entry becomes non-positive, factorization stops early and records the failure position.

Parameters:

get_elem – Function that returns matrix element at (i,j) position

Returns:

True if matrix is positive definite (all diagonal entries positive)

Return type:

bool

The factorization stores results in _storage: - Diagonal entries (D matrix) are stored in _storage[i,i] - Off-diagonal entries (L matrix) are stored in _storage[i,j] for j < i

Examples

>>> mat = np.array([[1.0, 0.5, 0.5], [0.5, 1.25, 0.75], [0.5, 0.75, 1.5]])
>>> ldl = LDLTMgr(3)
>>> ldl.factor(lambda i, j: mat[i, j])
True

factor_with_allow_semidefinite(get_elem: Callable[[int, int], float]) → bool

Performs LDLT factorization allowing for positive semi-definite matrices.

Similar to factor() but handles zero diagonal entries (indicating semi-definiteness) by restarting factorization from the next row.

Parameters:

get_elem – Function that returns matrix element at (i,j) position

Returns:

True if matrix is positive semi-definite (no negative diagonal entries)

Return type:

bool

This version is more tolerant of zero diagonal entries than factor(), which requires strictly positive entries for positive definiteness.

Examples

>>> mat = np.array([[1.0, 0.5, 0.5], [0.5, 1.25, 0.75], [0.5, 0.75, 1.5]])
>>> ldl = LDLTMgr(3)
>>> ldl.factor_with_allow_semidefinite(lambda i, j: mat[i, j])
True

factorize(mat: ndarray) → bool

Performs LDLT factorization on a given numpy array matrix.

This is a convenience wrapper around the factor() method that takes a matrix directly rather than an element-accessor function.

Parameters:

mat – A symmetric numpy array to be factorized

Returns:

True if matrix is positive definite, False otherwise

Return type:

bool

Examples

>>> mat = np.array([[1.0, 0.5, 0.5], [0.5, 1.25, 0.75], [0.5, 0.75, 1.5]])
>>> ldl = LDLTMgr(3)
>>> ldl.factorize(mat)
True

is_spd()

Checks if the matrix is symmetric positive definite (SPD).

Returns:

True if the matrix is SPD (pos[1] == 0), False otherwise

Return type:

bool

The check is based on whether any diagonal entry was non-positive during factorization, which would have set pos[1] to a non-zero value.

Examples

>>> mat = np.array([[1.0, 0.5, 0.5], [0.5, 1.25, 0.75], [0.5, 0.75, 1.5]])
>>> ldl = LDLTMgr(3)
>>> ldl.factorize(mat)
True
>>> ldl.is_spd()
True

pos: Tuple[int, int]

sqrt() → ndarray

Computes the upper triangular square root matrix R where A = R^T R.

Returns:

Upper triangular matrix R

Return type:

np.ndarray

Raises:

AssertionError – If matrix is not positive definite

This is essentially the Cholesky decomposition, computed from the LDLT factors without directly computing square roots until the final step.

Examples

>>> mat = np.array([[1.0, 0.5, 0.5], [0.5, 1.25, 0.75], [0.5, 0.75, 1.5]])
>>> ldl = LDLTMgr(3)
>>> ldl.factorize(mat)
True
>>> ldl.sqrt()
array([[1. , 0.5, 0.5],
       [0. , 1. , 0.5],
       [0. , 0. , 1. ]])

sym_quad(mat: ndarray)

Computes the quadratic form v^T M v using the witness vector.

Parameters:

mat – The matrix M to compute the quadratic form with

Returns:

The value of v^T M v where v is the witness vector

Return type:

float

Note: witness() must be called first to set up the witness vector. The computation uses only the submatrix where positive definiteness failed.

Examples

>>> mat = np.array([[1.0, 2.0, 3.0], [2.0, 3.5, 5.0], [3.0, 5.0, 6.0]])
>>> ldl = LDLTMgr(3)
>>> ldl.factorize(mat)
False
>>> ldl.pos
(0, 2)
>>> ldl.witness() # call this before sym_quad()
0.5
>>> ldl.wit
array([-2.,  1.,  0.])
>>> mat_b = np.array([[1.0, 0.5, 0.5], [0.5, 1.25, 0.75], [0.5, 0.75, 1.5]])
>>> ldl.sym_quad(mat_b)
3.25

wit: ndarray

witness() → float

Computes a witness vector proving the matrix is not positive definite.

Returns:

The negative eigenvalue (ep) showing v^T A v = -ep < 0

Return type:

float

Raises:

AssertionError – If called on a positive definite matrix

The witness vector is stored in self.wit and can be accessed after calling this method. The vector satisfies v^T A v < 0 for the failed submatrix.

Examples

>>> mat = np.array([[1.0, 2.0, 3.0], [2.0, 3.5, 5.0], [3.0, 5.0, 6.0]])
>>> ldl = LDLTMgr(3)
>>> ldl.factorize(mat)
False
>>> ldl.witness()
0.5

ellalgo.oracles.lmi0_oracle module

class ellalgo.oracles.lmi0_oracle.LMI0Oracle(mat_f)

Bases: object

Oracle for Linear Matrix Inequality (LMI) constraint: F(x) ⪰ 0

Solves the feasibility problem:

Find x ∈ ℝⁿ such that ∑_{k=1}^n F_k x_k ≽ 0

Where:

• F_k ∈ 𝕊^m (symmetric matrices) are given in mat_f

• x = [x_1, …, x_n]^T is the decision vector

• ≽ denotes positive semidefinite (PSD) constraint

The oracle uses LDLT factorization to verify PSD property and generates cutting planes for infeasible solutions.

assess_feas(x: ndarray) → Tuple[ndarray, float] | None

Assess feasibility of solution x against LMI constraint

Implementation Strategy: 1. Construct matrix F(x) = ∑ x_k F_k through element-wise computation 2. Attempt LDLT factorization:

• Success: F(x) is PSD (feasible) → return None

• Failure: F(x) not PSD → compute cutting plane (g, σ)

Parameters:

x (np.ndarray) – Candidate solution vector [x₁, …, xₙ]^T ∈ ℝⁿ

Returns:

• None if x is feasible (F(x) ≽ 0)

• Tuple (g, σ) representing cutting plane gᵀ(y - x) ≥ σ otherwise

Return type:

Optional[Cut]

ellalgo.oracles.lmi_old_oracle module

class ellalgo.oracles.lmi_old_oracle.LMIOldOracle(mat_f, mat_b)

Bases: OracleFeas

Oracle for Linear Matrix Inequality constraint.

This oracle solves the following feasibility problem:

find x s.t. (B − F * x) ⪰ 0

assess_feas(xc: ndarray) → Tuple[ndarray, float] | None

The assess_feas function assesses the feasibility of a given input array x and returns a Cut object if the feasibility is violated, otherwise it returns None.

Parameters:

x (np.ndarray) – An array of values that will be used in the calculation

Returns:

The function assess_feas returns an Optional[Cut].

ellalgo.oracles.lmi_oracle module

LMIOracle

This code defines a class called LMIOracle, which is designed to solve a specific type of mathematical problem known as a Linear Matrix Inequality (LMI) constraint. The purpose of this code is to determine if there exists a solution that satisfies the given constraint, and if not, to provide information about why it’s not feasible.

The LMIOracle class takes two inputs when it’s created: mat_f and mat_b. These are matrices that define the LMI constraint. mat_f is a list of numpy arrays, and mat_b is a single numpy array. These matrices represent the mathematical relationship that needs to be satisfied.

The main function in this class is called assess_feas, which takes a numpy array xc as input. This function tries to determine if xc is a feasible solution to the LMI constraint. If it is feasible, the function returns None. If it’s not feasible, it returns what’s called a “cut” - a tuple containing information about why the solution isn’t feasible.

To achieve its purpose, the code uses a technique called LDLT factorization. This is a way of breaking down a matrix into simpler parts, which can help determine if the matrix satisfies certain properties. The LDLTMgr class (which is used but not defined in this code snippet) handles this factorization.

The assess_feas function works by first creating a new matrix using the input xc and the original matrices mat_f and mat_b. It then tries to perform the LDLT factorization on this new matrix. If the factorization is successful, it means xc is a feasible solution, and the function returns None. If the factorization fails, it means xc is not feasible, and the function calculates and returns the “cut” information.

An important part of the logic is the get_elem function inside assess_feas. This function calculates each element of the new matrix based on the original matrices and the input xc. This is where the main mathematical operation of the LMI constraint is performed.

In summary, this code provides a way to check if a given solution satisfies a complex mathematical constraint, and if not, it provides information about why the solution doesn’t work. This could be useful in optimization problems where you’re trying to find a solution that satisfies certain mathematical conditions.

class ellalgo.oracles.lmi_oracle.LMIOracle(mat_f, mat_b)

Bases: OracleFeas

Oracle for Linear Matrix Inequality (LMI) constraint.

This oracle solves the following semidefinite feasibility problem:

find x

s.t. (B − ∑ F_i * x_i) ⪰ 0 [Matrix PSD constraint]

Where: - B is a constant symmetric matrix - F_i are coefficient matrices - x_i are decision variables - ⪰ 0 denotes positive semidefinite requirement

assess_feas(xc: ndarray) → Tuple[ndarray, float] | None

Assess feasibility of candidate solution xc against LMI constraint.

Implementation Steps: 1. Construct matrix M(xc) = B - ∑ F_i*xc_i 2. Perform LDLT factorization of M(xc) 3. If factorization succeeds (matrix PSD), return None (feasible) 4. If factorization fails (not PSD), return separating hyperplane (cut)

Parameters:

xc – Candidate solution vector xc ∈ ℝⁿ

Returns:

None if feasible, else (g, σ) where: - g ∈ ℝⁿ: subgradient of the violation - σ ∈ ℝ: measure of constraint violation

ellalgo.oracles.lowpass_oracle module

Lowpass Oracle

This code implements a Lowpass Oracle, which is used to design a low-pass filter for signal processing. A low-pass filter allows low-frequency signals to pass through while attenuating high-frequency signals. The main purpose of this code is to help optimize the design of such a filter by providing a way to assess whether a given set of filter coefficients meets certain specifications.

The code defines a class called LowpassOracle that takes several inputs when initialized:

1. ndim: The number of filter coefficients

2. wpass: The end of the passband (frequencies that should pass through)

3. wstop: The end of the stopband (frequencies that should be attenuated)

4. lp_sq: The lower bound for the squared magnitude response in the passband

5. up_sq: The upper bound for the squared magnitude response in the passband

6. sp_sq: The upper bound for the squared magnitude response in the stopband

The main outputs of this code are produced by two methods: assess_feas and assess_optim. These methods take a set of filter coefficients as input and determine whether they meet the specified requirements or how close they are to meeting them.

The LowpassOracle achieves its purpose through a series of checks on the frequency response of the filter. It uses a pre-computed spectrum matrix to efficiently calculate the frequency response at different points. The code then checks if the response falls within the specified bounds for the passband and stopband.

The important logic flow in this code involves iterating through different frequency points and checking the filter’s response at each point. If any violations of the specifications are found, the code returns information about the violation, which can be used to adjust the filter coefficients.

A key data transformation happening in this code is the conversion from filter coefficients to frequency response. This is done using the pre-computed spectrum matrix, which allows for efficient calculation of the response at many frequency points.

The code also includes a helper function called create_lowpass_case, which sets up a specific instance of the LowpassOracle with predefined parameters. This function can be used to quickly create a standard test case for filter design.

Overall, this code provides a tool for iteratively designing and optimizing low-pass filters by giving feedback on how well a set of coefficients meets the desired specifications. It’s part of a larger optimization process where the coefficients would be adjusted based on the feedback from this oracle until a satisfactory filter design is achieved.

class ellalgo.oracles.lowpass_oracle.LowpassOracle(ndim: int, wpass: float, wstop: float, lp_sq: float, up_sq: float, sp_sq: float)

Bases: OracleOptim

assess_feas(x: ndarray) → Tuple[ndarray, float | Tuple[float, float]] | None

Assess whether the given filter coefficients meet the design specifications.

This method checks the frequency response at various points in three bands: 1. Passband (0 to nwpass): Checks if response is within [lp_sq, up_sq] 2. Stopband (nwstop to end): Checks if response is below sp_sq and non-negative 3. Transition band (nwpass to nwstop): Checks if response is non-negative

Uses a round-robin approach to check different frequency points on each call to distribute the computational load across multiple iterations.

Parameters:

x (Arr) – The filter coefficients (autocorrelation coefficients)

Returns:

• None if all specifications are met

• A tuple containing:

  • The gradient of the violating constraint

  • The violation amount (or tuple of lower/upper violations)

Return type:

Optional[ParallelCut]

assess_optim(xc: ndarray, gamma: float)

Assess the optimality of the current filter coefficients for the stopband.

First checks feasibility using assess_feas. If feasible, returns information about the maximum response in the stopband which can be used to further optimize the filter design.

Parameters:

• xc (Arr) – The filter coefficients (autocorrelation coefficients)

• gamma (float) – The current best stopband attenuation value to beat

Returns:

A tuple containing:

• A tuple of (gradient, (lower, upper)) for the maximum stopband response

• The maximum stopband response value (or None if not feasible)

Return type:

tuple

idx1: int = 0

ellalgo.oracles.lowpass_oracle.create_lowpass_case(ndim=48)

Creates a standard low-pass filter design case with typical parameters.

Sets up a LowpassOracle instance with commonly used specifications: - Passband edge at 0.12π - Stopband edge at 0.20π - Passband ripple of ±0.025 dB - Stopband attenuation of 0.125

Parameters:

ndim (int, optional) – Number of filter coefficients. Defaults to 48.

Returns:

An initialized LowpassOracle instance with standard parameters

Return type:

LowpassOracle

ellalgo.oracles.profit_oracle module

Profit Oracle

This code defines several classes that implement oracles for profit maximization problems. An oracle, in this context, is a function that helps solve optimization problems by providing information about the feasibility and optimality of potential solutions.

The main class, ProfitOracle, is designed to solve a specific type of profit maximization problem. It takes as input parameters related to production (like price, scale, and limits) and output elasticities. The goal is to find the optimal input quantities that maximize profit, given certain constraints.

The ProfitOracle class has methods to assess the feasibility of a solution (assess_feas) and to find the optimal solution (assess_optim). These methods take as input a vector y (representing input quantities in log scale) and a gamma value (representing the current best solution). They output “cuts”, which are linear constraints that help narrow down the search for the optimal solution.

The code also includes two variations of the profit oracle:

1. ProfitRbOracle: This is a robust version of the profit oracle that can handle some uncertainty in the input parameters.

2. ProfitQOracle: This version deals with discrete (integer) input quantities, as opposed to continuous ones.

The main logic flow in these classes involves calculating various economic functions (like Cobb-Douglas production functions) and their gradients. The code uses these calculations to determine if a given solution is feasible and to guide the search towards the optimal solution.

The output of these oracles is typically a “cut” (a linear constraint) and possibly an updated best solution (gamma). These outputs are used by an external optimization algorithm (not shown in this code) to iteratively improve the solution until the optimal one is found.

For beginners, it’s important to understand that this code is implementing mathematical optimization techniques. While the details might be complex, the basic idea is to efficiently search for the best solution to a profit maximization problem, given certain constraints and economic relationships.

class ellalgo.oracles.profit_oracle.ProfitOracle(params: Tuple[float, float, float], elasticities: ndarray, price_out: ndarray)

Bases: OracleOptim

Oracle for a profit maximization problem using Cobb-Douglas production function.

This implementation follows the formulation from [Aliabadi and Salahi, 2013]:

Optimization problem:

max p(A y₁^α y₂^β) − v₁y₁ − v₂y₂ s.t. y₁ ≤ k

Where: - p(A y₁^α y₂^β): Cobb-Douglas production function in exponential form - p: Market price per unit - A: Production scale factor - α, β: Output elasticities (α + β ≤ 1 for constant returns to scale) - y: Input quantities (decision variables in log scale) - v: Input prices - k: Upper bound constraint for x₁

The oracle uses cutting plane methods to iteratively refine the solution space.

assess_feas(xc: ndarray, gamma: float) → Tuple[ndarray, float] | None

Feasibility assessment using round-robin constraint checking.

Implements: - Alternates between checking y₁ constraint (fn1) and optimality (fn2) - Returns first violated constraint found

Parameters:

• xc – Current solution point in log-space

• gamma – Current best profit estimate

Returns:

Cut (gradient, violation) if constraint violated None if all constraints satisfied

assess_optim(xc: ndarray, gamma: float) → Tuple[Tuple[ndarray, float], float | None]

Optimization assessment generating optimality cuts.

Workflow: 1. Check feasibility using assess_feas 2. If feasible, calculate new profit estimate and optimality cut 3. If infeasible, return feasibility cut

Parameters:

• y – Proposed solution point

• gamma – Current best profit value

Returns:

• Cut (gradient, beta)

• Updated gamma (None if infeasible)

Return type:

Tuple containing

elasticities: ndarray

fn1(x: ndarray, _: float) → float

Constraint function for y₁ ≤ k (in log-space).

Parameters:

x – Log-scale input vector [log(y₁), log(y₂)]

Returns:

x[0] - log(k) Positive values indicate constraint violation

Return type:

Constraint violation measure

fn2(x: ndarray, gamma: float) → float

Optimality condition function for profit maximization.

Computes: - Cobb-Douglas value in log-space: log(pA) + αlog(y₁) + βlog(y₂) - Variable costs: v₁y₁ + v₂y₂ - Optimality gap: log(γ + vy) - log_Cobb

Parameters:

• x – Log-scale input vector

• gamma – Current best profit estimate

Updates intermediate values used in gradient calculations

grad1(_: float) → ndarray

Gradient for y₁ ≤ k constraint.

Returns:

Gradient vector [1, 0] since ∂(x₀ - log_k)/∂x = (1, 0)

grad2(gamma: float) → ndarray

Gradient of optimality condition function.

Computes:

∇f = [v₁y₁/(γ+vy) - α, v₂y₂/(γ+vy) - β]

Parameters:

gamma – Current profit estimate used in denominator

Uses precomputed q (v₁y₁, v₂y₂) from last fn2 call

idx: int = -1

log_Cobb: float

log_k: float

log_pA: float

price_out: ndarray

q: ndarray

vy: float

class ellalgo.oracles.profit_oracle.ProfitQOracle(params, elasticities, price_out)

Bases: OracleOptimQ

Discrete profit oracle for integer input quantities.

Solves mixed-integer version of the profit maximization problem:

max p(A y₁^α y₂^β) − v₁y₁ − v₂y₂ s.t. x₁ ≤ k, x ∈ ℕ²

Uses continuous relaxation followed by rounding to nearest integer, with mechanisms to handle infeasible integer solutions.

assess_optim_q(xc: ndarray, gamma: float, retry: bool) → Tuple[Tuple[ndarray, float], ndarray, float | None, bool]

Discrete optimization assessment with rounding mechanism.

Workflow: 1. First try continuous solution (retry=False) 2. If infeasible, return feasibility cut 3. If feasible, round to nearest integer and check optimality 4. On retry (retry=True), check rounded solution optimality

Returns:

• Cut information

• Evaluation point (continuous or rounded)

• Updated gamma value

• Retry flag for integer feasibility check

Return type:

Tuple containing

xd: ndarray

class ellalgo.oracles.profit_oracle.ProfitRbOracle(params: Tuple[float, float, float], elasticities: ndarray, price_out: ndarray, vparams: Tuple[float, float, float, float, float])

Bases: OracleOptim

Robust profit oracle handling parameter uncertainty.

Implements robust optimization version from [Aliabadi and Salahi, 2013] considering uncertainties in: - Elasticity parameters (α, β) - Price parameters (p, v) - Production limit (k)

Uses interval-based uncertainty sets for robust constraint satisfaction.

assess_optim(xc: ndarray, gamma: float) → Tuple[Tuple[ndarray, float], float | None]

Robust optimization assessment accounting for parameter uncertainties.

Adjusts elasticities based on direction of uncertainty impact: - Decreases effective α, β when y > 0 (conservative adjustment) - Increases effective α, β when y ≤ 0 (worst-case scenario)