Large measurements b matrix


#1

Hi,
I would like to optimize a tikhonov objective function using a dense A matrix (m x n on the order of 1000 x 1000) on a series of measurements sampled over p time frames (b of size m x p , p on the order of 1000). A is ill-conditioned (cond# ~1E+12) as it comes from boundary element discretization of a Laplace operator.
I tried to pass the full time measurements b as a matrix but I got a different dimension size error. If I loop over time steps, the optimization takes very long. Is there a way to get the factorization of A once and keep re-using it to make the computation efficient/fast.
Thanks


(Mark L. Stone) #2

Can you show your correct (at least so you believe), but long-running code, as well your other code with error message it provided? Maybe include the output of whos for the code producing the dimension size. error message?


#3

Thanks Mark for responding. Any help is greatly appreciated. Here is the code I used.

%A;   %m x n matrix  
%PHI; %m x p matrix   
m=size(A,1);  
n=size(A,2);  
p=size(PHI,2);  
gamma = 0.04;  %regularization parameter fixed  
%%%%%% If I use full b matrix , error about size differences  
%%%%%% [[[??? Error using ==> cvx.plus at 45  Matrix dimensions must agree. 
%%%%%% Error in ==> cvx.minus at 21      
%%%%%% z = plus( x, y, true, cheat ); ]]] 

tic  
cvx_begin quiet
    b=PHI;
    variable x(n);
    minimize( norm(A*x-b) + gamma * norm(x,2) );  
cvx_end  
toc

%%%%%% If I use a loop, it takes ~ 3 seconds per loop iteration  
tic  
for k= 1:p
    cvx_begin quiet
        b=PHI(:,k);
        variable x(n);
        minimize( norm(A*x-b)+gamma * norm(x,2) );
    cvx_end   
    PHI_H(:,k)=x;  
end  
toc

(Mark L. Stone) #4

In your “full PHI matrix” code, Ax is m by 1, and b is my by p, so the dimensions don’t match to form Ax-b, hence MATLAB provides an error message. Your second code may be taking a long time because it incurs CVX overhead each time through the loop. Does each column of PHI correspond to a separate optimization problem, independent of the problems for other columns? If not, your second code is wrong. If the problems are independent, how about

cvx_begin
variable x(n,p)
expression Objective(p)
for k=1:p
Objective(k) = norm(A*x(:,k)-PHI(:,k)) + gamma * norm(x(:,k),2); 
end
minimize(sum(Objective))   
cvx_end

I have stacked p independent optimization problems on top of each other by summing the objectives, and I believe that the result will be that the solver will provide a solution to all at once, equivalent to if they had been presented as separate problems. So after CVX finishes, x(:,k) for k=1:p should contain the optimal x for the problem corresponding to the kth column of PHI.

This way, you should only incur the overhead of CVX once. I’ll leave any vectorization of the above as an exercise for you.

Edit: This edit serves only to bump the thread so that my comment below from Dec 22 will be evident.


Using CVX with multiple RHS
(Michael C. Grant) #5

In general, no, CVX provides no facility re-use calculations from iteration to iteration. CVX is designed for convenience, not for speed.


#6

Thanks Mark. I tried the stacking objective function which is a neat idea to go global optimization over all time steps however I got an out of memory error. I appreciate you help. The motivation behind using CVX is to use different objectives/constraints to see which would provide better results


(Mark L. Stone) #7

Then try stacking only some of the time steps together. Find a number of number of time steps you can stack without running out of memory (or requiring paging).