CVX long formulation time for large scale least square

[JeffWang-zxw20 (Jeff Wang)]

Hi,

I am using cvx to solve a least square problem. The objective function is:

minimise( norm(y - phi_u1*g1 - phi_u2*g2, 2) + gamma1*norm(g1,2) + gamma2*norm(g2,2) );

which is close to the standard L2-regularised least square problem. The size of y,u1,u2 are 2001x1, the size of phi_u1 and phi_u2 are 2001x400.

The problem is that the cvx does not produce any output for an hour. I guess cvx gets stuck at formulating the problem. Anyway to solve the problem? I know there is an analytical solution for this problem, but I plan to add LMI constraints to this objective function later, but cvx even can’t solve it without any constraints.

Any suggestion will be deeply appreciated!

Here is the full code

Out_ =  TwoFIR_cvx(u1_cvx,u2_cvx,y_cvx,para);

function Out_ = TwoFIR_cvx(u1,u2,y,para)
    %y = FIR(u1) + FIR(u2)
    if ~iscolumn(u1)
        error('need vec format input u1');
    end
    if ~iscolumn(u2)
        error('need vec format input u2');
    end
    if ~iscolumn(y)
        error('need vec format output y');
    end
    N = para.N;
    M1 = para.M1;
    M2 = para.M2;
    cvx_begin sdp 
        variables g1(M1,1) g2(M2,1);
        expressions phi_u1(N,M1) phi_u2(N,M2);  %can't declare here! otw bug!
        for i = 1:N
            for j = 1:M1
                if j>i
                    phi_u1(i,j) = 0.0;
                else
                    phi_u1(i,j) = u1(i-j+1);
                end
            end
        end

        for i = 1:N
            for j = 1:M2
                if j>i
                    phi_u2(i,j) = 0.0;
                else
                    phi_u2(i,j) = u2(i-j+1);
                end
            end
        end

        minimise( norm(y - phi_u1*g1 - phi_u2*g2, 2) + para.gamma1*norm(g1,2) + para.gamma2*norm(g2,2)  ); %minimise 2-norm


    cvx_end

    Out_.g1 = g1;
    Out_.G1 = filt(g1',1,para.Ts);
    Out_.phi_u1 = full(phi_u1);
    Out_.g2 = g2;
    Out_.G2 = filt(g2',1,para.Ts);
    Out_.phi_u2 = full(phi_u2);
    Out_.para = para;
end

[Mark_L_Stone]

There are two stages to getting the solution: Model formulation by CVX and model solving by the solver. Perhaps a lot of the time is due to CVX model formulation (during which there will be no output on the screen).

Vectorization should help tremendously with model formulation time.

How to vectorize most constraint loops in CVX

Using permute and repmat, most constraint loops in CVX can be eliminated. For example, the following: variables X(N) Y(M) Z(N); for n=1:N for m=1:M if n~=m && n^3>m^2 X(n)>=Y(m)+Z(n); % (1) end end end The first step is to expand. (1) actually means there are M*N constraints, so X is expanded to repmat(X,1,M), Y is expanded to repmat(Y,1,N), Z(n) is expanded to repmat(Z,1,M), according to their original dimensions. Now, repmat(X,1,M) is a N-by-M matrix, rep…

[JeffWang-zxw20 (Jeff Wang)]

Hi Mark,

Thank you so much for your response! I will have a try of your suggestion!

So you think this has nothing to do with the solver choice ( I actually tried the other free solver but the same thing happens)? Also, switching to cvxpy or cvx Julia will not help at all (planning to do so)?

Thanks

[Mark_L_Stone]

I have no idea what will happen if the solver ever gets called. Vectorization should also help in CVXPY or Julia.

Perhaps you can try a small version of the problem and see what happens.