"Mark " <markcouwenberg@atgmaildot.com> wrote in message <ganstp$1je$1@fred.mathworks.com>...
> Dear All,
>
> I want to find an optimal solution for a linear set of equations with more unknowns than equations, ultimately up to 3 equations and some 410 variables.
> Lets start simple with 2 equations, 3 unknowns
> coefficient matrix:
> A=[1 1 0; 0 0 1];
> b = [10;10];
>
> Solution represent a set of forces which have a certain upper limit, e.g.: ub=[20 20 20];
>
> I do the following:
> 1. find particular solution:
> Fp = A\b;
> 2. find null solution:
> F0 = null(A);
>
> Now all possible solutions are given by:
> F = Fp + F0 * const
>
> and the target is to find a constant which gives an optimal solution. Optimal I define as follows:
> We are talking about force equilibrium, F represent forces which have a certain maximum (ub) and minimum (lb).
> My solution is optimal if the sum of the square of the ratios F/Fmax is minimal:
> obj = sum( (F./ub).^2 )
> In this simple example it is easy to handcalculate that this is obtained with F = [5 5 10]
>
> I use optimisation tool fmincon. For the F0 multiplier "const" I add a new variable to my solution and add a column of zeros to A which gives Aeq.
> Using the following piece of code:
> >>>
> %%%Add zeros after A
> Aeq(:,numF+1:numX) = 0;
> beq = b;
> x0=[0; 0; 0; 0];
> options = optimset('LargeScale','off');
> [x,fval,exitflag,output] = fmincon(@(x) Fobjfun(x, Fp, F0, ub),x0,[],[],Aeq,beq,lb,ub,[],options)
> >>>
> With Fobjfun:
> >>>function obj = Fobjfun(x, Fp, F0, ub)
> x(1:3) = Fp(i) + F0(i)*x(4);
> obj = sum( (x(1:numF)./ub(1:numF)).^2 );
> >>>
>
> This works, it gives the answer: x=[5;5;10;7.0711] and obj 0.375.
> Next:
> I lower one constraint:
> ub = [4; 20; 20]
>
> Now I get the answer:
> >>>
> x =
>
> 3.9255
> 6.0745
> 10.0000
> 13.5982
> obj =
>
> 0.4904
> >>>
> Verification by simply trying a range of possible C's (x4) shows indeed that x4 shall be 13.5982 and the corresponding obj value is 0.4904. However the corresponding forces are not gives in x(1:3). I have to add a final line which reconstructs forces as is done in the objective funtion. This seems strange to me. x4 does not correspond with x1:x3. If I reconstruct forces with F=Fp + F0 * x4 than a set of forces is given which is differnt from x1x3 but gives the correct (and indeed optimal) objective value.
> Questions:
> 1. Anyone understands what is going wrong here? Or is it not going wrong?
> 2. Anyone has better ideas in finding the optimal solution from a set of linear equations with more unknowns than equations?
>
> Any help is greatly appreciated,
>
> Mark
>
>
Dear All,
I found the solution. The reason I use a multiplier for my null solution is the way I previously "solved" these kind of problems: simply by trying a huge range of possibilties of multipliers. However using decent optimisation technology this is unnecessary. I can enter my coefficient matrix directly as Aeq. There is no need for adding an extra variable for every null space vector. So now my code would be:
>>>
A = [1 1 0 0 0;...
0 0 1 1 1;...
1 1 0.5 1 1];
beq = [10;10;10];
ub = [20; 20; 20; 20; 20];
lb = [20; 20; 20; 20; 20];
Aeq = A;
x0=[0; 0; 0; 0; 0];
options = optimset('LargeScale','off');
[x,fval,exitflag,output] = fmincon(@(x) Fobjfun(x, ub),x0,[],[],Aeq,beq,lb,ub,[],options)
>>>
and the objective function:
>>>
function obj = Fobjfun(x, ub)
obj = sum( (x./ub).^2 );
>>>
This works. However I am still curious about what happens in the previous posted method. Any ideas?
Brgds,
Mark
