MATLAB: Most efficient way to sum anti-diagonal elements

diagonalsindexing

Inside an ODE that I am solving (ode15s), I need to do the following. Let the state vector by y (Nx1), and for some fixed, sparse symmetric A (NxN), I need to have that
 y' = [sum of anti-diagonal elements of (diag(y)*A*diag(y))] + f(t)
for some forcing f. The problem is, for large N (10k + ) this is pretty slow as the solver takes many steps. Currently I have calculated a logical index mask (outside the ode) that gives me the elements I want, and my code (inside the d.e.) is:
A = spdiags(y,0,N,N)*A*spdiags(y,0,N,N); %overwrite A in-place with y*A*y
working_matrix=sparse([],[],[],2*N-1,N,N*ceil(N/2)); % sparse allocation
working_matrix(index_map)=A; %strip antidiagonals to columns
this gives me working_matrix, which has the antidiagonal elements of (diag(y)*A*diag(y) which I can just sum over. However, 99% of my runtime is spent on the line "working_matrix(index_map)=A". Any speedup on this line would save me a lot of time. "index_map" is a (2*N-1)xN logical array that pulls out the correct elements, based on this work here.
Is there a better way? I can pull of the antidiagonal elements of A outside the solver and pass the matrix that has the antidiagonal elements as rows, but then I still need the same construction applied to y*y' to get the matching elements of y – unless there is a better way to do this?
I am running R2015b if that matters.

Best Answer

"all of the elements on any NE-SW diagonal"
a = randi(50,6);
[m,n] = size(a);
idx = hankel(1:m,m:(n-1)+m);
out = accumarray(idx(:),a(:));
Related Question