Optimization with Lagrangian multiplier

Question

I have a function, fd, which I hope to optimize using an equality constraint. Here is the function:

fd = 224 * d1 + 84 * d2 + d1 * d2 - 2 * d1^2 - d2^2

and constraint:

3*d1 + d2 = 280

Someone has solved this problem by incorporating the constraint into the function using a Lagrangian multiplier. The new function Ld.lambda is:

Ld.lambda = 224 * d1 + 84 * d2 + d1 * d2 - 2 * d1^2 - d2^2 + lambda * (280 - 3*d1 - d2)

They used the following R code to obtain the values of d1, d2 and lambda that maximize, I think, the new function, Ld.lambda:

# Lagrangian multiplier method
obj_c <- function(param)
                 {x1     <- param[1]
                  x2     <- param[2]
                  lambda <- param[3]
                  L      <- 224*x1+84*x2+x1*x2-2*x1^2-x2^2+lambda*(280-3*x1-x2)
                  L
         }

# gradient function or partial derivatives of Ld.lambda
grad_c <- function(dec)
                  {x1     <- dec[1]
                   x2     <- dec[2]
                   lambda <- dec[3]
                   g1     <- 224+x2-4*x1-3*lambda
                   g2     <- 84+x1-2*x2-lambda
                   g3     <- 280-3*x1-x2
                   c(g1,g2,g3)
           }

# linear solver to estimate optimality conditions
a <- matrix(data<-c(4,-1,3,-1,2,1,3,1,0),nrow=3,ncol=3,byrow=T)
b <- c(224,84,280)
dec <- solve(a,b)
dec
# [1] 69 73  7

The optimal values are estimated to be:

d1 = 69
d2 = 73
lambda = 7

However, if I attempt a quasi-exhaustive search in R the optimal values appear to be:

d1 = 69
d2 = 73
lambda = -infinity to +infinity

Apparently because:

(280 - 3*d1 - d2) = 0 when d1 = 69 and d2 = 73, in which case the value of lambda is irrelevant.

Here is the R code for the quasi-exhaustive search:

my.data <- expand.grid(x1 = seq(0, 200, 1), x2 = seq(0, 200, 1), x3 = seq(0, 200, 1))
head(my.data)
dim(my.data)

d1     <- my.data[,1]
d2     <- my.data[,2]
lambda <- my.data[,3]

Fd <- 224 * d1 + 84 * d2 + d1 * d2 - 2 * d1^2 - d2^2 + lambda * (280 - 3*d1 - d2)

new.data <- data.frame(Fd = Fd, d1 = d1, d2 = d2, lambda = lambda)
head(new.data)

# Impose constraint
new.data <- new.data[(3 * new.data$d1 + new.data$d2) == 280, ]

# identify values of d1, d2 and lambda that maximize Fd with the constraint
head(new.data[new.data$Fd == max(new.data$Fd),], 10)
tail(new.data[new.data$Fd == max(new.data$Fd),], 10)

Why does the Lagrangian multiplier method estimate the optimal value of lambda = 7 and the quasi-exhaustive search estimate the optimal value of lambda = -infinity to +infinity?

What am I doing incorrectly or misunderstanding?

Answer 1

Bottom line up front: You should be solving a system of equations involving the Lagrange multiplier, not maximizing it. Details are below.

Details: The Lagrange multiplier method consists of solving the Karush-Kuhn-Tucker (KKT) conditions for first order optimum, i.e., actually, stationary point. See https://en.wikipedia.org/wiki/Karush%E2%80%93Kuhn%E2%80%93Tucker_conditions . Without looking at 2nd order conditions, you can not determine whether the solution is actually an optimum. In your case, because you have no inequality constraints, the KKT conditions take the simpler form (special case) shown in https://en.wikipedia.org/wiki/Lagrange_multiplier, and are the same for minimization and maximization; so in addition to the possibility of finding a saddle point, you could be finding a local maximum or local minimum.

Note that although your problem falls under the simpler special case described in the 2nd link, the first link provides a much better theoretical explanation as to what is going on, which is needed to understand what is going on in the 2nd link.

However, because your objective function is concave (and sufficiently differentiable) and your constraint is linear, thereby satisfying regularity conditions (constraint qualification), KKT conditions are necessary and sufficient for the solution to be a global maximum. Indeed, your objective function is strictly concave, so in this case, the KKT conditions have a single solution, which is the global maximum.

So I think the first method is solving the KKT conditions, which in your case is (just) a system of linear equations.

If you had any inequality constraints, you would need to use the more general form of KKT conditions shown in the first link, and those aren't necessarily easy to solve. In practice, many numerical optimizers for constrained optimization set up subproblems to drive toward a solution to the KKT conditions, which in "typical" cases do not admit a closed form or non-iterative solution, although they do in some cases. See Variance of weighted mean greater than unweighted mean for an example in which KKT conditions involving inequality constraints can be solved in closed form, although I didn't show the details of how i did it.