Fast matrix copying and array appending

Hi,

I would like to know if there exists a fast way to copy a subpart of a matrix into a new matrix and to append an array to the same matrix. For instance, something like this: if true A=[A(2:end,2:end),B;B',b]; end where A is nxn matrix, B is nx1 array and b is a scalar.

Thank you Best Regards Salvo

Answer 1

[Jan]

I'm convinced that this is a good idea already:

A = [A(2:end, 2:end), B; B', b];

You save some micro-seconds by avoiding end, because the explicit definition if the size is faster:

[s1, s2] = size(A);
A = [A(2:s1, 2:s2), B; B', b];

In my measurements this is remarkably slower:

A(1:s1-1, 1:s2-1) = A(2:s1, 2:s2);
A(s1,     1:s2-1) = B';
A(1:s1-1, s2)     = B;
A(s1, s2)         = b;

This might save some temporary memory, but the indexing needs more time.

[EDITED]

void mexFunction(int nlhs, mxArray *plhs[], int nrhs, const mxArray *prhs[]) {
  double *A, *B, *Bp, b, *R;
  mwSize i, n, n_m1;

    % Get inputs: 
    A = mxGetPr(prhs[0]);
    B = mxGetPr(prhs[1]);
    b = mxGetScalar(prhs[2]);
    n = mxGetM(prhs[0]);

    % Create outputs:
    plhs[0] = mxCreateDoubleMatrix(n, n, mxReal);
    R = mxGetPr(plhs[0]);

    Bp = B;
    n_m1 = n - 1;
    for (i = 0; i < n_m1; i++) {
      memcpy(R, A+1, n_m1 * sizeof(double));
      R[n_m1] = *Bp++;
      R += n;
      A += n;
    }

    memcpy(R, B, n_m1 * sizeof(double));
    R[n_m1] = b;
  }

Use this with great care and after exhaustive testing only! This is written in the forum's editor and neither compiled nor tested. Do you know how to compile a C-Mex file?

[nedo nodo]

So I cannot improve the speed of my code, is it right?

[José-Luis]

It depends whether you want your data in a single matrix. Depending on the application, the size of your matrix, your access patterns, etc., it might be more efficient to use a cell array. It is very difficult to know what will be faster without actually testing it. And even then, the results will only be true for your computer.

[Jan]

@nedo nodo: I do not believe, that this assignment needs the most time in your program. Optimizing a not relevant part of the code does not lead to a significant acceleration of the whole program. So please use the profiler to find the bottlenecks, if you did not do this already.

You could create a C-MEX function for the creation of the matrix, but it is not sure if this is faster at all, because most of Matlab's indexing functions are very efficient already. If you find out, that this operation needs 90% of the overall processing time, I could offer to create a small C-code.

[nedo nodo]

@Jan Simone I have already done the profiling of my code and for this reason I have said that it is fundamental to improve this row. Your solution without the esplicit use of end is sligthly faster. Thank you

Answer 2

[Matt J]

If this really is the performance critical part of the code, you should consider redesigning so that A(2:end,2:end) gets prepended instead of appended. It leads to orders of magnitude speed difference, as shown below. The redesign should be trivial if you are doing likely operations like A*x. Just do A*[x(end);x(1:end-1)] instead.

    N=3000;
    A=rand(N);
    B=rand(N-1,1);
    b=1;

#########

    tic %Appending
     A=[A(2:end,2:end),B;B',b];
    toc
    %Elapsed time is 0.086694 seconds.

############

    tic%Prepending

      A(1)=b;
      A(2:end,1)=B;
      A(1,2:end)=B';
    toc
    %Elapsed time is 0.000218 seconds.

[Matt J]

You could also use my ProdCascade class so that even though you're prepending A, it will automatically pre- and post-permute things you want to multiply with A to give the same I/O behavior as if you had appended to A. The class introduces some overhead, but it's still close to three times as fast as your original approach when doing A*x operations:

tic %Appending
 A=[A(2:end,2:end),B;B',b];
 z1=A*x;
toc

Elapsed time is 0.079278 seconds.

P=circshift(speye(N),[1,0]);
AA=ProdCascade({P.',A,P});

tic;

    A(1)=b;
    A(2:end,1)=B;
    A(1,2:end)=B';
    AA{2}=A;
    z2=AA*x;
  toc;

Elapsed time is 0.027279 seconds.