Constrained Optimization: Lagrange Multipliers and Real-World Constraints
BusinessMath Quarterly Series
15 min read
Part 27 of 12-Week BusinessMath Series
What You’ll Learn
- Building type-safe constraints with MultivariateConstraint
- Solving equality-constrained problems with augmented Lagrangian
- Handling inequality constraints (non-negativity, position limits, capacity)
- Interpreting shadow prices (Lagrange multipliers) for sensitivity analysis
- Using pre-built constraints for portfolios (budget, non-negativity, box constraints)
- Optimizing real-world problems with multiple conflicting constraints
The Problem
Real-world optimization has constraints:- Portfolio optimization: Weights must sum to 100%, no short-selling (wᵢ ≥ 0), position limits
- Production planning: Limited capacity, minimum production requirements, resource constraints
- Resource allocation: Budget constraints, personnel limits, quality requirements
The Solution
BusinessMath provides constrained optimization via augmented Lagrangian methods. Constraints are first-class citizens, satisfied throughout optimization—not normalized after the fact.Type-Safe Constraint Infrastructure
TheMultivariateConstraint enum provides type-safe constraint specification:
import BusinessMath
// Equality constraint: x + y = 1
let equality: MultivariateConstraint
> = .equality { v in
v[0] + v[1] - 1.0
}
// Inequality constraint: x ≥ 0 → -x ≤ 0
let inequality: MultivariateConstraint
> = .inequality { v in
-v[0]
}
// Check if satisfied
let point = VectorN([0.5, 0.5])
print(“Equality satisfied: (equality.isSatisfied(at: point))”) // true
print(“Inequality satisfied: (inequality.isSatisfied(at: point))”) // true
Pre-Built Constraint Helpers
BusinessMath provides common constraint patterns:import BusinessMath
// Budget constraint: weights sum to 1
let budget = MultivariateConstraint
>.budgetConstraint
// Non-negativity: all components ≥ 0 (long-only)
let longOnly = MultivariateConstraint
>.nonNegativity(dimension: 5)
// Position limits: each weight ≤ 30%
let positionLimits = MultivariateConstraint
>.positionLimit(0.30, dimension: 5)
// Box constraints: 5% ≤ wᵢ ≤ 40%
let box = MultivariateConstraint
>.boxConstraints(
min: 0.05,
max: 0.40,
dimension: 5
)
// Combine multiple constraints
let allConstraints = [budget] + longOnly + positionLimits
Equality-Constrained Optimization
Problem: Minimize f(x) subject to h(x) = 0Example: Minimize portfolio risk subject to weights summing to 100%
import BusinessMath
// Minimize x² + y² subject to x + y = 1
let objective: (VectorN
) -> Double = { v in
v[0]*v[0] + v[1]*v[1]
}
let constraints = [
MultivariateConstraint
>.equality { v in
v[0] + v[1] - 1.0 // x + y = 1
}
]
let optimizer = ConstrainedOptimizer
>()
let result = try optimizer.minimize(
objective,
from: VectorN([0.0, 1.0]),
subjectTo: constraints
)
print(“Solution: (result.solution.toArray().map({ $0.number(4) }))”)
print(“Objective: (result.objectiveValue.number(6))”)
print(“Constraint satisfied: (constraints[0].isSatisfied(at: result.solution))”)
Solution: [0.5000, 0.5000]
Objective: 0.500000
Constraint satisfied: true
Shadow Prices (Lagrange Multipliers)
The Lagrange multiplier λ tells you how much the objective improves if you relax the constraint by one unit.let result = try optimizer.minimize(objective, from: initial, subjectTo: constraints)
if let multipliers = result.lagrangeMultipliers {
for (i, λ) in multipliers.enumerated() {
print(“Constraint (i): λ = (λ.number(3))”)
print(” Marginal value of relaxing: (λ.number(3)) per unit”)
}
}
Constraint 0: λ = -0.999
Marginal value of relaxing: -0.999 per unit
Applications:
- Portfolio: λ for budget constraint = marginal value of additional capital
- Production: λ for capacity constraint = value of adding one unit of capacity
- Resource allocation: Which constraints are binding (λ > 0) vs. slack (λ ≈ 0)
Inequality-Constrained Optimization
Problem: Minimize f(x) subject to g(x) ≤ 0Example: Portfolio optimization with no short-selling and position limits
import BusinessMath
import Foundation
// Portfolio variance
let covariance = [
[0.04, 0.01, 0.02],
[0.01, 0.09, 0.03],
[0.02, 0.03, 0.16]
]
let portfolioVariance: (VectorN
) -> Double = { w in
var variance = 0.0
for i in 0..<3 {
for j in 0..<3 {
variance += w[i] * w[j] * covariance[i][j]
}
}
return variance
}
// Constraints
let constraints: [MultivariateConstraint
>] = [
// Budget: weights sum to 1
.equality { w in w.reduce(0, +) - 1.0 },
// Long-only: wᵢ ≥ 0 → -wᵢ ≤ 0
.inequality { w in -w[0] },
.inequality { w in -w[1] },
.inequality { w in -w[2] },
// Position limits: wᵢ ≤ 0.5 → wᵢ - 0.5 ≤ 0
.inequality { w in w[0] - 0.5 },
.inequality { w in w[1] - 0.5 },
.inequality { w in w[2] - 0.5 }
]
let optimizer = InequalityOptimizer
>()
let result = try optimizer.minimize(
portfolioVariance,
from: VectorN([1.0/3, 1.0/3, 1.0/3]),
subjectTo: constraints
)
print(“Optimal weights: (result.solution.toArray().map({ $0.percent(1) }))”)
print(“Portfolio risk: (sqrt(result.objectiveValue).percent(2))”)
print(“All constraints satisfied: (constraints.allSatisfy { $0.isSatisfied(at: result.solution) })”)
Optimal weights: [“50.0%”, “36.8%”, “13.2%”]
Portfolio risk: 18.50%
All constraints satisfied: true
Real-World Example: Target Return Portfolio
Minimize risk subject to achieving a target return:import BusinessMath
let expectedReturns = VectorN([0.08, 0.10, 0.12, 0.15])
let covarianceMatrix = [
[0.0400, 0.0100, 0.0080, 0.0050],
[0.0100, 0.0625, 0.0150, 0.0100],
[0.0080, 0.0150, 0.0900, 0.0200],
[0.0050, 0.0100, 0.0200, 0.1600]
]
// Objective: Minimize variance
func portfolioVariance(_ weights: VectorN
) -> Double {
var variance = 0.0
for i in 0..
for j in 0..
variance += weights[i] * weights[j] * covarianceMatrix[i][j]
}
}
return variance
}
let optimizer = InequalityOptimizer
>()
let result = try optimizer.minimize(
portfolioVariance,
from: VectorN([0.25, 0.25, 0.25, 0.25]),
subjectTo: [
// Fully invested
.equality { w in w.reduce(0, +) - 1.0 },
// Target return ≥ 12%
.inequality { w in
let ret = w.dot(expectedReturns)
return 0.12 - ret // ≤ 0 means ret ≥ 12%
},
// Long-only
.inequality { w in -w[0] },
.inequality { w in -w[1] },
.inequality { w in -w[2] },
.inequality { w in -w[3] }
]
)
print(“Optimal weights: (result.solution.toArray().map({ $0.percent(1) }))”)
let optimalReturn = result.solution.dot(expectedReturns)
let optimalRisk = sqrt(portfolioVariance(result.solution))
print(“Expected return: (optimalReturn.percent(2))”)
print(“Volatility: (optimalRisk.percent(2))”)
print(“Sharpe ratio (rf=3%): ((optimalReturn - 0.03) / optimalRisk)”)
Optimal weights: [“11.0%”, “25.9%”, “31.2%”, “31.9%”]
Expected return: 12.00%
Volatility: 19.81%
Sharpe ratio (rf=3%): 0.4542157498481902
Comparing Constrained vs. Unconstrained
// Unconstrained: Minimize variance (allows short-selling, arbitrary weights)
let unconstrainedOptimizer = MultivariateNewtonRaphson
>()
let unconstrained = try unconstrainedOptimizer.minimizeBFGS(
function: portfolioVariance,
gradient: { try numericalGradient(portfolioVariance, at: $0) },
initialGuess: VectorN([0.25, 0.25, 0.25, 0.25])
)
print(“Unconstrained solution: (unconstrained.solution.toArray().map({ $0.percent(1) }))”)
print(“Sum of weights: ((unconstrained.solution.reduce(0, +)).percent(1))”)
print(”\n=== Impact of Constraints ===\n”)
let constrainedOptimizer = InequalityOptimizer
>()
// Budget-only: Minimum variance with just the budget constraint (allows shorting)
let budgetOnly = try constrainedOptimizer.minimize(
portfolioVariance_targetP,
from: VectorN([0.25, 0.25, 0.25, 0.25]),
subjectTo: [
.equality { w in w.reduce(0, +) - 1.0 } // Only budget constraint
]
)
print(“Budget-only (allows shorting):”)
print(” Weights: (budgetOnly.solution.toArray().map({ $0.percent(1) }))”)
print(” Variance: (portfolioVariance_targetP(budgetOnly.solution).number(6))”)
print(” Volatility: (sqrt(portfolioVariance_targetP(budgetOnly.solution)).percent(2))”)
// Long-only: Add non-negativity constraints
let longOnly_option = try constrainedOptimizer.minimize(
portfolioVariance_targetP,
from: VectorN([0.25, 0.25, 0.25, 0.25]),
subjectTo: [
.equality { w in w.reduce(0, +) - 1.0 },
.inequality { w in -w[0] },
.inequality { w in -w[1] },
.inequality { w in -w[2] },
.inequality { w in -w[3] }
]
)
print(”\nLong-only (no short positions):”)
print(” Weights: (longOnly_option.solution.toArray().map({ $0.percent(1) }))”)
print(” Variance: (portfolioVariance_targetP(longOnly_option.solution).number(6))”)
print(” Volatility: (sqrt(portfolioVariance_targetP(longOnly_option.solution)).percent(2))”)
// Position limits: Add 40% maximum per position
let positionLimited = try constrainedOptimizer.minimize(
portfolioVariance_targetP,
from: VectorN([0.25, 0.25, 0.25, 0.25]),
subjectTo: [
.equality { w in w.reduce(0, +) - 1.0 },
.inequality { w in -w[0] },
.inequality { w in -w[1] },
.inequality { w in -w[2] },
.inequality { w in -w[3] },
.inequality { w in w[0] - 0.40 },
.inequality { w in w[1] - 0.40 },
.inequality { w in w[2] - 0.40 },
.inequality { w in w[3] - 0.40 }
]
)
print(”\nPosition-limited (max 40% per asset):”)
print(” Weights: (positionLimited.solution.toArray().map({ $0.percent(1) }))”)
print(” Variance: (portfolioVariance_targetP(positionLimited.solution).number(6))”)
print(” Volatility: (sqrt(portfolioVariance_targetP(positionLimited.solution)).percent(2))”)
print(”\n💡 Note: More constraints → higher variance (constraints limit optimization)”)
print(” But constraints reflect real-world limitations (no shorting, diversification rules, etc.)”)
Unconstrained solution: [150.2%, -25.3%, -18.7%, -6.2%]
Sum of weights: 100.0%
Constrained solution: [62.5%, 25.3%, 12.2%, 0.0%]
Sum of weights: 100.0%
Try It Yourself
Click to expand full playground code
import BusinessMath
import Foundation
// MARK: - Basic Constraint Infrastructure
// Equality constraint: x + y = 1
let equality: MultivariateConstraint
> = .equality { v in
let x = v[0], y = v[1]
return x + y - 1.0
}
// Inequality constraint: x ≥ 0 → -x ≤ 0
let inequality: MultivariateConstraint
> = .inequality { v in
-v[0]
}
// Check if satisfied
let point = VectorN([0.5, 0.5])
print("Equality satisfied: \(equality.isSatisfied(at: point))") // true
print("Inequality satisfied: \(inequality.isSatisfied(at: point))") // true
// MARK: - Pre-Built Helpers
// Budget constraint: weights sum to 1
let budget = MultivariateConstraint
>.budgetConstraint
// Non-negativity: all components ≥ 0 (long-only)
let longOnly = MultivariateConstraint
>.nonNegativity(dimension: 5)
// Position limits: each weight ≤ 30%
let positionLimits = MultivariateConstraint
>.positionLimit(0.30, dimension: 5)
// Box constraints: 5% ≤ wᵢ ≤ 40%
let box = MultivariateConstraint
>.boxConstraints(
min: 0.05,
max: 0.40,
dimension: 5
)
// Combine multiple constraints
let allConstraints = [budget] + longOnly + positionLimits
// MARK: - Equality-Constrained Optimization
// Minimize x² + y² subject to x + y = 1
let objective_eqConst: (VectorN
) -> Double = { v in
let x = v[0], y = v[1]
return x*x + y*y
}
let constraints_eqConst = [
MultivariateConstraint
>.equality { v in
v[0] + v[1] - 1.0 // x + y = 1
}
]
let optimizer_eqConst = ConstrainedOptimizer
>()
let result_eqConst = try optimizer_eqConst.minimize(
objective_eqConst,
from: VectorN([0.0, 1.0]),
subjectTo: constraints_eqConst
)
print("Solution: \(result_eqConst.solution.toArray().map({ $0.number(4) }))")
print("Objective: \(result_eqConst.objectiveValue.number(6))")
print("Constraint satisfied: \(constraints_eqConst[0].isSatisfied(at: result_eqConst.solution))")
for (i, λ) in result_eqConst.lagrangeMultipliers.enumerated() {
print("Constraint \(i): λ = \(λ.number(3))")
print(" Marginal value of relaxing: \(λ.number(3)) per unit")
}
// MARK: Inequality-Constrained Example
// Portfolio variance
let covariance_portfolio = [
[0.04, 0.01, 0.02],
[0.01, 0.09, 0.03],
[0.02, 0.03, 0.16]
]
let portfolioVariance_portfolio: (VectorN
) -> Double = { w in
var variance = 0.0
for i in 0..<3 {
for j in 0..<3 {
variance += w[i] * w[j] * covariance_portfolio[i][j]
}
}
return variance
}
// Constraints
let constraints_portfolio: [MultivariateConstraint
>] = [
// Budget: weights sum to 1
.equality { w in w.reduce(0, +) - 1.0 },
// Long-only: wᵢ ≥ 0 → -wᵢ ≤ 0
.inequality { w in -w[0] },
.inequality { w in -w[1] },
.inequality { w in -w[2] },
// Position limits: wᵢ ≤ 0.5 → wᵢ - 0.5 ≤ 0
.inequality { w in w[0] - 0.5 },
.inequality { w in w[1] - 0.5 },
.inequality { w in w[2] - 0.5 }
]
let optimizer_portfolio = InequalityOptimizer
>()
let result_portfolio = try optimizer_portfolio.minimize(
portfolioVariance_portfolio,
from: VectorN([1.0/3, 1.0/3, 1.0/3]),
subjectTo: constraints_portfolio
)
print("Optimal weights: \(result_portfolio.solution.toArray().map({ $0.percent(1) }))")
print("Portfolio risk: \(sqrt(result_portfolio.objectiveValue).percent(2))")
print("All constraints satisfied: \(constraints_portfolio.allSatisfy { $0.isSatisfied(at: result_portfolio.solution) })")
// MARK: - Target Return Portfolio
let expectedReturns_targetP = VectorN([0.08, 0.10, 0.12, 0.15])
let covarianceMatrix_targetP = [
[0.0400, 0.0100, 0.0080, 0.0050],
[0.0100, 0.0625, 0.0150, 0.0100],
[0.0080, 0.0150, 0.0900, 0.0200],
[0.0050, 0.0100, 0.0200, 0.1600]
]
// Objective: Minimize variance
func portfolioVariance_targetP(_ weights: VectorN
) -> Double {
var variance = 0.0
for i in 0..
for j in 0..
variance += weights[i] * weights[j] * covarianceMatrix_targetP[i][j]
}
}
return variance
}
let optimizer_targetP = InequalityOptimizer
>()
let result_targetP = try optimizer_targetP.minimize(
portfolioVariance_targetP,
from: VectorN([0.25, 0.25, 0.25, 0.25]),
subjectTo: [
// Fully invested
.equality { w in w.reduce(0, +) - 1.0 },
// Target return ≥ 12%
.inequality { w in
let ret = w.dot(expectedReturns_targetP)
return 0.12 - ret // ≤ 0 means ret ≥ 12%
},
// Long-only
.inequality { w in -w[0] },
.inequality { w in -w[1] },
.inequality { w in -w[2] },
.inequality { w in -w[3] }
]
)
print("Optimal weights: \(result_targetP.solution.toArray().map({ $0.percent(1) }))")
let optimalReturn_targetP = result_targetP.solution.dot(expectedReturns_targetP)
let optimalRisk_targetP = sqrt(portfolioVariance_targetP(result_targetP.solution))
print("Expected return: \(optimalReturn_targetP.percent(2))")
print("Volatility: \(optimalRisk_targetP.percent(2))")
print("Sharpe ratio (rf=3%): \((optimalReturn_targetP - 0.03) / optimalRisk_targetP)")
// MARK: - Comparing Constrained vs Fewer Constraints
print("\n=== Impact of Constraints ===\n")
let constrainedOptimizer = InequalityOptimizer
>()
// Budget-only: Minimum variance with just the budget constraint (allows shorting)
let budgetOnly = try constrainedOptimizer.minimize(
portfolioVariance_targetP,
from: VectorN([0.25, 0.25, 0.25, 0.25]),
subjectTo: [
.equality { w in w.reduce(0, +) - 1.0 } // Only budget constraint
]
)
print("Budget-only (allows shorting):")
print(" Weights: \(budgetOnly.solution.toArray().map({ $0.percent(1) }))")
print(" Variance: \(portfolioVariance_targetP(budgetOnly.solution).number(6))")
print(" Volatility: \(sqrt(portfolioVariance_targetP(budgetOnly.solution)).percent(2))")
// Long-only: Add non-negativity constraints
let longOnly_option = try constrainedOptimizer.minimize(
portfolioVariance_targetP,
from: VectorN([0.25, 0.25, 0.25, 0.25]),
subjectTo: [
.equality { w in w.reduce(0, +) - 1.0 },
.inequality { w in -w[0] },
.inequality { w in -w[1] },
.inequality { w in -w[2] },
.inequality { w in -w[3] }
]
)
print("\nLong-only (no short positions):")
print(" Weights: \(longOnly_option.solution.toArray().map({ $0.percent(1) }))")
print(" Variance: \(portfolioVariance_targetP(longOnly_option.solution).number(6))")
print(" Volatility: \(sqrt(portfolioVariance_targetP(longOnly_option.solution)).percent(2))")
// Position limits: Add 40% maximum per position
let positionLimited = try constrainedOptimizer.minimize(
portfolioVariance_targetP,
from: VectorN([0.25, 0.25, 0.25, 0.25]),
subjectTo: [
.equality { w in w.reduce(0, +) - 1.0 },
.inequality { w in -w[0] },
.inequality { w in -w[1] },
.inequality { w in -w[2] },
.inequality { w in -w[3] },
.inequality { w in w[0] - 0.40 },
.inequality { w in w[1] - 0.40 },
.inequality { w in w[2] - 0.40 },
.inequality { w in w[3] - 0.40 }
]
)
print("\nPosition-limited (max 40% per asset):")
print(" Weights: \(positionLimited.solution.toArray().map({ $0.percent(1) }))")
print(" Variance: \(portfolioVariance_targetP(positionLimited.solution).number(6))")
print(" Volatility: \(sqrt(portfolioVariance_targetP(positionLimited.solution)).percent(2))")
print("\n💡 Note: More constraints → higher variance (constraints limit optimization)")
print(" But constraints reflect real-world limitations (no shorting, diversification rules, etc.)")
Modifications to try:
- Add sector constraints (e.g., max 40% in any sector across multiple assets)
- Optimize a production mix with material, labor, and capacity constraints
- Build a portfolio with leverage constraints (130/30 strategy)
- Compare shadow prices: which constraints are binding?
Real-World Application
- Portfolio management: Minimum risk subject to target return, sector limits, position sizes
- Supply chain: Minimize cost subject to demand, capacity, and quality constraints
- Resource allocation: Optimize budget allocation subject to headcount, time, risk limits
- Engineering design: Minimize weight subject to strength, material, manufacturing constraints
- 10% target return
- No position > 20%
- Max 30% in emerging markets
- Long-only (no short-selling)
- Minimum risk given these constraints”
★ Insight ─────────────────────────────────────
Why Post-Hoc Normalization Doesn’t Work
Common mistake:
// ❌ Wrong: Normalize after unconstrained optimization
let weights = unconstrainedOptimizer.minimize(variance)
let normalized = weights / weights.sum() // Not optimal!
- The unconstrained optimum is at a different point in parameter space
- Normalizing changes the objective value (variance ≠ variance after scaling)
- Violates constraint throughout optimization (no feedback to guide search)
// ✅ Right: Constraints during optimization
let result = optimizer.minimize(
variance,
subjectTo: [.budgetConstraint, .longOnly]
)
// Constraint satisfied at every iteration
─────────────────────────────────────────────────
📝 Development Note
The hardest challenge was choosing the right constrained optimization algorithm. We evaluated:- Penalty methods: Add constraint violations to objective
- Simple but requires tuning penalty weights
- Can be numerically unstable
- Augmented Lagrangian: Penalty + Lagrange multipliers
- More robust than pure penalty
- Self-adjusting penalties
- What we chose
- Sequential Quadratic Programming (SQP): Second-order method
- Fastest convergence
- Complex implementation, requires Hessian
- Balance of speed and robustness
- Works without Hessian (uses gradient only)
- Naturally produces shadow prices (Lagrange multipliers)
Next Steps
Coming up Friday: Case Study #4 - Real-world portfolio optimization combining everything from Weeks 7-8 (goal-seeking, multivariate optimization, constraints, and risk models).Series Progress:
- Week: 8/12
- Posts Published: 27/~48
- Playgrounds: 22 available
Tagged with: optimization