Implementing matrix multiplication constraints in fmincon

I am attempting to solve a chemistry problem which involves optimizing the concentrations of 5 chemicals at 32 points along a reaction. For a given point along the reaction, the weighted concentrations of the 5 chemicals must sum to a known scalar, where the weighting involves multiplying each of the concentrations by an appropriate scalar (i.e., multiply the first concentration by 1, the second by 1, the third by 2, etc.). The problem can be formulated in matrix multiplication as [1,1,2,3,4]*C = b, where C is the 5x32 concentration matrix and b is the vector representing the sum of the concentrations for each point along the reaction.

I would like to use fmincon to optimize C for an objective function I have (which takes C as a matrix) subject to the constraint outlined above. However, I keep on getting errors saying that Aeq (the [1,1,2,3,4] weighting vector) must have 160 elements. I understand from other posts that linear constraints appear to be enforced after splitting in the matrix into vectors, but my problem does not allow me to optimize each point along the reaction individually. Is there any way a matrix multiplication constraint can be successfully implement in fmincon?

Thanks,

Nathanael

Answer 1

[John D'Errico]

You DO have 160 unknowns from what you have said. The 5 concentrations, at each of 32 positions. It is just that you are thinking of them as a 5x32 array.

Apparently you have a linear equality constraint at EACH of the 32 locations, that links the concentrations at each location.

So build a linear constraint matrix and right hand side. (Thus the A and b inputs to fmincon.) A will be a 32 by 160 matrix, that is mainly zeros. (You can build it using blkdiag.) b, the right hand side vector here, will be a 32x1 vector, composed of the required weighted sum at each location.

So, perhaps something like this:

C = repmat({[1 1 2 3 4]},1,32);
A = blkdiag(C{:});
spy(A)

Answer 2

[Alan Weiss]

Your control variables C is a 5x32 matrix, which represents 160 elements. As the documentation describes, each row of a linear constraint matrix must be the same length as the number of control variables.

I assume that you want your linear constraint to be an equality constraint, of the form

Aeq*x = beq

where I write the conventional x instead of your notation C. In this form, each row of Aeq multiplies the vector x, and this vector is the column vector

x = C(:);

So how do you turn your problem into the MATLAB form? You need to represent Aeq in such as way that the result Aeq*x is a column vector of 32 entries that should equal the vector b. I suggest using a sparse matrix representation as follows:

Aeq = spalloc(32,160,160);
bb =  [1,1,2,3,4];
for ii = 1:32
    Aeq(ii,(5*ii-4:5*ii)) = bb;
end

Then use a 32-element vector for beq and everything should work.

Good luck,

Alan Weiss

MATLAB mathematical toolbox documentation

[Alan Weiss]

I just noticed that John answered this before I did. Actually, I think we started to answer at the same time, but I got interrupted so finished formulating my answer after his was posted, but I didn't see it. You will see that the constraint matrices we suggest are identical.

Alan Weiss

MATLAB mathematical toolbox documentation

[John D'Errico]

This happens a lot. :) Hey, at least we agreed!

[Nathanael Kazmierczak]

Thank you both very much for your detailed help!