Discover MakerZone

MATLAB and Simulink resources for Arduino, LEGO, and Raspberry Pi

Learn more

Discover what MATLAB® can do for your career.

Opportunities for recent engineering grads.

Apply Today

Thread Subject:
Optimising set of linear equations

Subject: Optimising set of linear equations

From: Mark

Date: 16 Sep, 2008 09:06:01

Message: 1 of 2

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 4-10 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 hand-calculate 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 x1-x3 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

Subject: Optimising set of linear equations

From: Mark

Date: 16 Sep, 2008 13:24:02

Message: 2 of 2

"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 4-10 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 hand-calculate 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 x1-x3 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

Tags for this Thread

What are tags?

A tag is like a keyword or category label associated with each thread. Tags make it easier for you to find threads of interest.

Anyone can tag a thread. Tags are public and visible to everyone.

Contact us