Bilinear Product Linearization

Guide for linearizing bilinear products (x × y) in optimization models.

Overview

Bilinear products appear frequently in optimization:

  • Revenue: price × quantity

  • Energy: power × time

  • Logistics: distance × flow

  • Finance: return × allocation

The LXBilinearLinearizer class automatically selects the appropriate linearization technique based on variable types.

Variable Type Combinations

The linearization method depends on the types of the two variables:

Type 1

Type 2

Method

Complexity

Binary

Binary

AND logic

3 constraints

Binary

Continuous

Big-M

4 constraints

Binary

Integer

Big-M

4 constraints

Continuous

Continuous

McCormick

4 constraints

Continuous

Integer

McCormick

4 constraints

Integer

Integer

McCormick

4 constraints

Binary × Binary (AND Logic)

Mathematical Formulation

For z = x × y where x, y ∈ {0, 1}:

\[\begin{split}z \leq x \\ z \leq y \\ z \geq x + y - 1 \\ z \in \{0, 1\}\end{split}\]

Logical Interpretation:

  • z = 1 if and only if both x = 1 and y = 1

  • Implements logical AND operation

Example

from lumix import LXModel, LXVariable
from lumix.linearization.techniques import LXBilinearLinearizer
from lumix.linearization.config import LXLinearizerConfig
from lumix.nonlinear.terms import LXBilinearTerm

# Define binary variables
is_open = LXVariable[Facility, int]("is_open").binary()
is_served = LXVariable[Customer, int]("is_served").binary()

# Create bilinear term: facility_open AND customer_served
term = LXBilinearTerm(
    var1=is_open,
    var2=is_served,
    coefficient=1.0
)

# Linearize
config = LXLinearizerConfig()
linearizer = LXBilinearLinearizer(config)
z = linearizer.linearize_bilinear(term)

# z represents is_open AND is_served

Use Cases

  • Facility Location: Open facility AND serve customer

  • Scheduling: Assign shift AND assign driver

  • Logic: Boolean AND operations in constraints

Binary × Continuous (Big-M Method)

Mathematical Formulation

For z = b × x where b ∈ {0, 1}, x ∈ [L, U]:

\[\begin{split}z \leq M \cdot b \\ z \geq m \cdot b \\ z \leq x - m(1 - b) \\ z \geq x - M(1 - b) \\ z \in [m, M]\end{split}\]

where M = upper bound of x, m = lower bound of x.

Behavior:

  • When b = 0: z = 0 (forced by first two constraints)

  • When b = 1: z = x (forced by last two constraints)

Example

# Binary variable: is facility open?
is_open = (
    LXVariable[Facility, int]("is_open")
    .binary()
    .indexed_by(lambda f: f.id)
    .from_data(facilities)
)

# Continuous variable: flow amount
flow = (
    LXVariable[Route, float]("flow")
    .continuous()
    .bounds(lower=0, upper=1000)
    .indexed_by(lambda r: r.id)
    .from_data(routes)
)

# Product: flow only if facility is open
# z = is_open × flow
term = LXBilinearTerm(var1=is_open, var2=flow, coefficient=1.0)

config = LXLinearizerConfig(big_m_value=1000)
linearizer = LXBilinearLinearizer(config)
z = linearizer.linearize_bilinear(term)

Choosing Big-M Value

Critical: M must be large enough but not too large.

Guidelines:

  1. Use Variable Bounds:

    flow = LXVariable[Route, float]("flow").bounds(lower=0, upper=1000)
# M = 1000 (upper bound)

  2. Problem Knowledge:

    # For normalized probabilities [0, 1]
config = LXLinearizerConfig(big_m_value=1.0)

# For prices up to $10,000
config = LXLinearizerConfig(big_m_value=10000)

  3. Conservative Multiplier:

    max_flow = 1000
config = LXLinearizerConfig(big_m_value=10 * max_flow)

Validation:

# After solving, check that auxiliary variable doesn't hit M
solution = optimizer.solve(linearized_model)

for var_name, value in solution.variables.items():
    if "aux_prod" in var_name:
        assert value < 0.99 * config.big_m_value, \
            f"Variable {var_name} hit Big-M bound!"

Use Cases

  • Conditional Flow: Flow only if route is open

  • Pricing: Revenue only if product is sold

  • Scheduling: Work hours only if worker is assigned

Continuous × Continuous (McCormick Envelopes)

Mathematical Formulation

For z = x × y where x ∈ [xL, xU], y ∈ [yL, yU]:

\[\begin{split}z \geq x_L \cdot y + y_L \cdot x - x_L \cdot y_L \\ z \geq x_U \cdot y + y_U \cdot x - x_U \cdot y_U \\ z \leq x_L \cdot y + y_U \cdot x - x_L \cdot y_U \\ z \leq x_U \cdot y + y_L \cdot x - x_U \cdot y_L\end{split}\]

Geometric Interpretation:

These four constraints form the convex hull of the bilinear function over the box [xL, xU] × [yL, yU].

Note

McCormick envelopes provide the tightest possible linear relaxation of the bilinear term over the given bounds.

Example

# Continuous variable: price ($/unit)
price = (
    LXVariable[Product, float]("price")
    .continuous()
    .bounds(lower=10, upper=100)
    .indexed_by(lambda p: p.id)
    .from_data(products)
)

# Continuous variable: quantity (units)
quantity = (
    LXVariable[Product, float]("quantity")
    .continuous()
    .bounds(lower=0, upper=1000)
    .indexed_by(lambda p: p.id)
    .from_data(products)
)

# Product: revenue = price × quantity
term = LXBilinearTerm(var1=price, var2=quantity, coefficient=1.0)

config = LXLinearizerConfig(
    mccormick_tighten_bounds=True  # Improve relaxation
)
linearizer = LXBilinearLinearizer(config)
revenue = linearizer.linearize_bilinear(term)

Importance of Tight Bounds

Tighter bounds = stronger relaxation = faster solving

Example of bound impact:

# Wide bounds (weak relaxation)
x = LXVariable[Model, float]("x").bounds(lower=0, upper=1000)
y = LXVariable[Model, float]("y").bounds(lower=0, upper=1000)
# McCormick creates large convex hull

# Tight bounds (strong relaxation)
x = LXVariable[Model, float]("x").bounds(lower=10, upper=50)
y = LXVariable[Model, float]("y").bounds(lower=5, upper=20)
# McCormick creates tight convex hull

Bound Tightening:

Enable automatic bound tightening:

config = LXLinearizerConfig(
    mccormick_tighten_bounds=True,  # Apply preprocessing
    auto_detect_bounds=True  # Use variable bounds
)

Use Cases

  • Revenue: price × quantity

  • Energy: power × time

  • Economics: elasticity × demand

  • Engineering: force × displacement

Advanced: Multi-dimensional Products

For products of more than two variables, apply linearization recursively:

# For x × y × z:
# Step 1: w = x × y (linearize)
# Step 2: result = w × z (linearize)

from lumix.linearization.techniques import LXBilinearLinearizer

config = LXLinearizerConfig()
linearizer = LXBilinearLinearizer(config)

# First product: w = x × y
term_xy = LXBilinearTerm(var1=x, var2=y, coefficient=1.0)
w = linearizer.linearize_bilinear(term_xy)

# Second product: result = w × z
term_wz = LXBilinearTerm(var1=w, var2=z, coefficient=1.0)
result = linearizer.linearize_bilinear(term_wz)

# Add all auxiliary elements to model
for var in linearizer.auxiliary_vars:
    model.add_variable(var)
for constraint in linearizer.auxiliary_constraints:
    model.add_constraint(constraint)

Complete Example

Production Planning with Revenue Maximization

from lumix import LXModel, LXVariable, LXConstraint
from lumix.linearization import LXLinearizer, LXLinearizerConfig
from lumix.solvers import LXOptimizer

# Define data
@dataclass
class Product:
    id: str
    min_quantity: float
    max_quantity: float
    min_price: float
    max_price: float

products = [
    Product("A", 0, 1000, 10, 100),
    Product("B", 0, 500, 20, 150),
]

# Define variables
quantity = (
    LXVariable[Product, float]("quantity")
    .continuous()
    .indexed_by(lambda p: p.id)
    .bounds_func(lambda p: (p.min_quantity, p.max_quantity))
    .from_data(products)
)

price = (
    LXVariable[Product, float]("price")
    .continuous()
    .indexed_by(lambda p: p.id)
    .bounds_func(lambda p: (p.min_price, p.max_price))
    .from_data(products)
)

# Build model
model = LXModel("production")

# Objective: maximize revenue (bilinear: price × quantity)
# This will be automatically linearized
revenue_expr = price * quantity  # Nonlinear!
model.maximize(revenue_expr)

# Constraints
# Total quantity constraint
model.add_constraint(
    LXConstraint("total_quantity")
    .expression(quantity)
    .le()
    .rhs(1200)
)

# Configure linearization
config = LXLinearizerConfig(
    default_method=LXLinearizationMethod.MCCORMICK,
    mccormick_tighten_bounds=True,
    verbose_logging=True
)

# Solve with linearization
optimizer = LXOptimizer().use_solver("glpk")
solver_capability = optimizer.get_capability()

linearizer = LXLinearizer(model, solver_capability, config)
linearized_model = linearizer.linearize_model()

# Get statistics
stats = linearizer.get_statistics()
print(f"Linearization added:")
print(f"  - {stats['auxiliary_variables']} variables")
print(f"  - {stats['auxiliary_constraints']} constraints")

# Solve
solution = optimizer.solve(linearized_model)
print(f"Optimal revenue: ${solution.objective_value:,.2f}")

Performance Tips

  1. Always Provide Bounds

    # Good: McCormick will be tight
x = LXVariable[Model, float]("x").bounds(lower=0, upper=100)

# Bad: Wide default bounds, weak relaxation
x = LXVariable[Model, float]("x")  # No bounds!

  2. Enable Bound Tightening

    config = LXLinearizerConfig(mccormick_tighten_bounds=True)

  3. Use Appropriate Big-M

    # Problem-specific, not default
config = LXLinearizerConfig(big_m_value=1000)  # Not 1e6!

  4. Prefer Native Quadratic When Available

    # For Gurobi/CPLEX with quadratic support,
# solver may handle bilinear natively without linearization

See Also