Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Do variables need to be greater than zero? #3

Open
mengchaoheng opened this issue Mar 24, 2024 · 1 comment
Open

Do variables need to be greater than zero? #3

mengchaoheng opened this issue Mar 24, 2024 · 1 comment

Comments

@mengchaoheng
Copy link

Usually, the standard form of linear programming requires variables to be greater than zero. Does this constraint need to be integrated into the A matrix?

@mengchaoheng
Copy link
Author

I have a LP problem in matlab code is:

% If we use the Standard Forms for Linear Programming Problems
% min c'x subj. to A*x <= b
%                    0 <= x
%% so we have to reformula the direction-preserving control allocation
% problem to:
% min z=[0; -1]'[u; a]   s.t.  [B -v][u; a] = 0
%                                umin <= u <= umax
%                                   0 <= a
% and set x=u-umin, then
% min z=[0; -1]'[x; a]   s.t.  [B -v][x; a] = -B*umin
%                                        x <= umax-umin
%                                        0 <= x 
%                                        0 <= a
% min z=[0; -1]'[x; a]   s.t.  [B -v][x; a] <= -B*umin
%                              [-B v][x; a] <=  B*umin
%                              [I  0][x; a] <= umax-umin
%                                        -x <= 0 
%                                        -a <= 0
% set X=[x;a]
l1=0.148;l2=0.069;k=3;
B=k*[-l1     0       l1     0;
     0      -l1     0       l1;
     l2    l2    l2    l2];
umin=ones(m,1)*(-20)*pi/180;
umax=ones(m,1)*20*pi/180;
[k,m] = size(B);
v=[0.2;0.1;0.1];
A=[B -v;-B v; eye(m) zeros(m,1)];
b=[-B*umin; B*umin; umax-umin];
c=[zeros(m,1); -1];
[X,fval,exitflag,output,lambda]= linprog(c,A,b,[],[],zeros(2*m+1,1),[]);
u = X(1:m)+umin;
a = X(m+1);
% Scale down u if a>1
if a>1
    u = u/a;
end
u

I have sole this by matlab function linprog.
and then solve this by SDLP, for SDLP , A, b and c is

% A=[B -v;-B v; eye(m) zeros(m,1);-eye(m) zeros(m,1);zeros(1,m) -1]
% b=[-B*umin; B*umin; umax-umin;0;0;0;0;0]
% c=[zeros(m,1); -1]
% or
A=[B -v;-B v; eye(m) zeros(m,1)]
b=[-B*umin; B*umin; umax-umin]
c=[zeros(m,1); -1]

setup A, b and c to sdlp_example and then run :

build git:(main) ✗ make          
[ 50%] Building CXX object CMakeFiles/sdlp_example.dir/example/sdlp_example.cpp.o
[100%] Linking CXX executable sdlp_example
[100%] Built target sdlp_examplebuild git:(main) ✗ ./sdlp_example
optimal sol: 0.155327 0.426713   0.6981   0.6981  1.20496
u: -0.193773 0.0776133     0.349     0.349
optimal obj: -1.20496

get the solution:

u_SDLP=[-0.193773; 0.0776133;     0.349;     0.349];

the optimal obj is the same, but B*u_SDLP != v, u != u_SDLP. why? @ZhepeiWang

the sdlp_example code is

#include <iostream>

#include "sdlp/sdlp.hpp"

using namespace std;
using namespace Eigen;

int main(int argc, char **argv)
{
    // int m = 2 * 7;
    int m = 10;
    Eigen::Matrix<double, 5, 1> x;        // decision variables
    Eigen::Matrix<double, 5, 1> c;        // objective coefficients
    Eigen::Matrix<double, -1, 5> A(m, 5); // constraint matrix
    Eigen::VectorXd b(m);                 // constraint bound
    A <<   -0.4440,         0,              0.4440,    0,        -0.2000,
            0,             -0.4440,         0,         0.4440,   -0.1000,
            0.2070,         0.2070,         0.2070,    0.2070,   -0.1000,
            0.4440,         0,             -0.4440,    0,         0.2000,
            0,              0.4440,         0,        -0.4440,    0.1000,
           -0.2070,        -0.2070,        -0.2070,   -0.2070,    0.1000,
            1.0000,         0,              0,         0,         0,
            0,              1.0000,         0,         0,         0,
            0,              0,              1.0000,    0,         0,
            0,              0,              0,         1.0000,    0;
    b <<0, 0, 0.2890, 0, 0, -0.2890, 0.6981, 0.6981, 0.6981, 0.6981;
    c << 0.0, 0.0, 0.0, 0.0, -1.0;
    double minobj = sdlp::linprog<5>(c, A, b, x);
    Eigen::Matrix<double, 4, 1> u_={x(0)-0.3491, x(1)-0.3491, x(2)-0.3491, x(3)-0.3491};
    Eigen::Matrix<double, 4, 1> u=u_;
    if(minobj>=1)
    {
        for(int i=0;i<4;i++)
        {
            u(i)=u_(i)/minobj;
        }
    }
    std::cout << "optimal sol: " << x.transpose() << std::endl;
    std::cout << "u: " << u.transpose() << std::endl;
    std::cout << "optimal obj: " << minobj << std::endl;
    return 0;
}

Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
Labels
None yet
Projects
None yet
Development

No branches or pull requests

1 participant