Thread Subject: Replacing the diagonal of a square matrix

Hey guys,

Is there a way to replace the diagonal of a square matrix with something else without a loop?

Like if we have,

a=[1 2 3;4 5 6;7 8 9]

Replace the diagonal, so [1 5 9] with Inf, so we get,

a=[Inf 2 3;4 Inf 6;7 8 Inf]

Thanks!

"Diego Zegarra"
> a=[1 2 3;4 5 6;7 8 9]
> Replace the diagonal, so [1 5 9] with Inf, so we get,
> a=[Inf 2 3;4 Inf 6;7 8 Inf]

one of the solutions

     a=[
          1,2,3
          4,5,6
          7,8,9
     ]
     a(eye(size(a))~=0)=inf
%{
   Inf 2 3
     4 Inf 6
     7 8 Inf
%}

us

For square matrix a:

a(linspace(1,numel(a),length(a))) = Inf

or

a(1:length(a)+1:numel(a)) = Inf








_fin(pob\gn!efi_cyo4_G[n_yy]c[h?ijyya[[iy[[jbbi]!Shi[\syg:y

"Diego Zegarra" <diegozbb@gmail.com> wrote in message <gopkna$9cf$1@fred.mathworks.com>...
> Hey guys,
>
> Is there a way to replace the diagonal of a square matrix with something else without a loop?
>
> Like if we have,
>
> a=[1 2 3;4 5 6;7 8 9]
>
> Replace the diagonal, so [1 5 9] with Inf, so we get,
>
> a=[Inf 2 3;4 Inf 6;7 8 Inf]
>
> Thanks!

a=[1 2 3;4 5 6;7 8 9];
n = size(a,1);

a(1:(n+1):end) = inf
a =
   Inf 2 3
     4 Inf 6
     7 8 Inf

HTH,
John

"Diego Zegarra" <diegozbb@gmail.com> wrote in message <gopkna$9cf$1@fred.mathworks.com>...
> Hey guys,
>
> Is there a way to replace the diagonal of a square matrix with something else without a loop?
----------------------------

I don't think there's any reason to. I think a for-loop will be the fastest way.

See for example, the following comparison wiht Matt Fig's approach

a=rand(3000); v=linspace(1,numel(a),length(a));

>> tic;a(v)=inf;toc
Elapsed time is 0.000356 seconds.

>> tic;for ii=1:length(a), a(ii)=inf; end,toc
Elapsed time is 0.000143 seconds.

"Matt " <xys@whatever.com> wrote in message <gopml5$mmp$1@fred.mathworks.com>...

> I don't think there's any reason to. I think a for-loop will be the fastest way.
>
> See for example, the following comparison wiht Matt Fig's approach
>
> a=rand(3000); v=linspace(1,numel(a),length(a));
>
> >> tic;a(v)=inf;toc
> Elapsed time is 0.000356 seconds.
>
> >> tic;for ii=1:length(a), a(ii)=inf; end,toc
> Elapsed time is 0.000143 seconds.


Your For loop did not replace the elements on the diagonal. A fairer comparison is the folowing:

a=rand(3000);
b = a;
b(2) = 5; % Force the copy.
tic;b(1:length(b)+1:end)=inf;toc
b(2) = a(2);
tic;for ii=1:length(a)+1:numel(a), a(ii)=inf; end,toc

all(a(:)==b(:)) % True;
Elapsed time is 0.000221 seconds.
Elapsed time is 0.000180 seconds.

So indeed you are correct. The For loop is slightly faster.

"Matt " <xys@whatever.com> wrote in message <gopml5$mmp$1@fred.mathworks.com>...
> "Diego Zegarra" <diegozbb@gmail.com> wrote in message <gopkna$9cf$1@fred.mathworks.com>...
> > Hey guys,
> >
> > Is there a way to replace the diagonal of a square matrix with something else without a loop?
> ----------------------------
>
> I don't think there's any reason to. I think a for-loop will be the fastest way.
>
> See for example, the following comparison wiht Matt Fig's approach
>
> a=rand(3000); v=linspace(1,numel(a),length(a));
>
> >> tic;a(v)=inf;toc
> Elapsed time is 0.000356 seconds.
>
> >> tic;for ii=1:length(a), a(ii)=inf; end,toc
> Elapsed time is 0.000143 seconds.

By the way, the loop above is NOT the same as
replacing the diagonal elements!!!!!!!!! See that
a(ii) does not touch those diagonal elements.
Instead, this replaced the first column of a. Even
that was far better done without a loop.

But seriously, do you honestly think so a for loop
will even be remotely competitive? I'll try it. First
of all why in god's name use linspace, when colon
is perfect? These are integer indices after all. The
calls to linspace and numel just create extra function
call overhead. They are not necessary.

a=rand(3000);

tic;a(1:(n+1):end)=inf;toc
Elapsed time is 0.000736 seconds.

tic,for i = 1:n,a(i,i)= inf;end,toc
Elapsed time is 0.020839 seconds.

It looks to me as if the for loop was approximately
a factor of 30 times slower than the direct index.

John

"Matt Fig" <spamanon@yahoo.com> wrote in message <gopnac$89p$1@fred.mathworks.com>...
> "Matt " <xys@whatever.com> wrote in message <gopml5$mmp$1@fred.mathworks.com>...
>
> > I don't think there's any reason to. I think a for-loop will be the fastest way.
> >
> > See for example, the following comparison wiht Matt Fig's approach
> >
> > a=rand(3000); v=linspace(1,numel(a),length(a));
> >
> > >> tic;a(v)=inf;toc
> > Elapsed time is 0.000356 seconds.
> >
> > >> tic;for ii=1:length(a), a(ii)=inf; end,toc
> > Elapsed time is 0.000143 seconds.
>
>
> Your For loop did not replace the elements on the diagonal. A fairer comparison is the folowing:
>
> a=rand(3000);
> b = a;
> b(2) = 5; % Force the copy.
> tic;b(1:length(b)+1:end)=inf;toc
> b(2) = a(2);
> tic;for ii=1:length(a)+1:numel(a), a(ii)=inf; end,toc
>
> all(a(:)==b(:)) % True;
> Elapsed time is 0.000221 seconds.
> Elapsed time is 0.000180 seconds.
>
> So indeed you are correct. The For loop is slightly faster.

This one may save time, due to the use of single
element indexing, as opposed to dual subscripts.

But it very much depends upon the matlab release.
It is possible they may have made some JIT
enhancements, because even the single index
loop is no faster in my release.

tic,for i = 1:n,a(i,i)= inf;end,toc
Elapsed time is 0.021719 seconds.

tic,for i = 1:(n+1):n^2,a(i)= inf;end,toc
Elapsed time is 0.021161 seconds.

tic;a(1:(n+1):end)=inf;toc
Elapsed time is 0.000712 seconds.

john

"John D'Errico" <woodchips@rochester.rr.com> wrote in message
> It looks to me as if the for loop was approximately
> a factor of 30 times slower than the direct index.
>
> John


It looks like you don't have the accelerator turned on. When I ran your code I got the same results as those I posted earlier. The For loop IS faster. Not that it is prettier or even faster by much, but the numbers speak for themselves.

r2007a WinVista 32

Just before I posted this I checked something that may be of interest. At the command line, I get results similar to what you get. In a script file, I get the for loops being faster. When I make the script a function (which allows Matlab's in-place optimizations) the For loop is MUCH faster. Instead of relative timings of 221/180, I get 222/148.

Interesting....

"Matt Fig" <spamanon@yahoo.com> wrote in message <goppit$c2j$1@fred.mathworks.com>...
> "John D'Errico" <woodchips@rochester.rr.com> wrote in message
> > It looks to me as if the for loop was approximately
> > a factor of 30 times slower than the direct index.
> >
> > John
>
>
> It looks like you don't have the accelerator turned on. When I ran your code I got the same results as those I posted earlier. The For loop IS faster. Not that it is prettier or even faster by much, but the numbers speak for themselves.
>
> r2007a WinVista 32
>
> Just before I posted this I checked something that may be of interest. At the command line, I get results similar to what you get. In a script file, I get the for loops being faster. When I make the script a function (which allows Matlab's in-place optimizations) the For loop is MUCH faster. Instead of relative timings of 221/180, I get 222/148.
>
> Interesting....

FYI, the JIT does not operate on code typed in at the command line.

James Tursa

Well, I'm finding that both for loops (I fixed my earlier mistake) and while loops are slightly slower than colon notation, but not by much. I do have the accelerator turned on, I think.

a=rand(3000); N=length(a); M=numel(a);


%colon notation
tic; a(1:N+1:M) = Inf; toc
Elapsed time is 0.000241 seconds.

%while loop
 tic; ii=1; while ii<=N, a(ii,ii)=inf; ii=ii+1; end,toc
Elapsed time is 0.000269 seconds.

%for loop
tic; for ii=1:N, a(ii,ii)=inf; end,toc
Elapsed time is 0.000269 seconds.

"James Tursa" <aclassyguywithaknotac@hotmail.com> wrote in message <gopupq$841$1@fred.mathworks.com>...
> "Matt Fig" <spamanon@yahoo.com> wrote in message <goppit$c2j$1@fred.mathworks.com>...
> > "John D'Errico" <woodchips@rochester.rr.com> wrote in message
> > > It looks to me as if the for loop was approximately
> > > a factor of 30 times slower than the direct index.
> > >
> > > John
> >
> >
> > It looks like you don't have the accelerator turned on. When I ran your code I got the same results as those I posted earlier. The For loop IS faster. Not that it is prettier or even faster by much, but the numbers speak for themselves.
> >
> > r2007a WinVista 32
> >
> > Just before I posted this I checked something that may be of interest. At the command line, I get results similar to what you get. In a script file, I get the for loops being faster. When I make the script a function (which allows Matlab's in-place optimizations) the For loop is MUCH faster. Instead of relative timings of 221/180, I get 222/148.
> >
> > Interesting....
>
> FYI, the JIT does not operate on code typed in at the command line.
>
> James Tursa


In light of this, I revised my tests from my last post by sticking them in a function.

Now the while loop is the fastest, as I would expect since it does not force the creation of an index array.

%%%%%%%%%%%%%%%%%%%%%%
function diagtest

a=rand(3000); N=length(a); M=numel(a);

tic; a(1:N+1:M) = Inf; toc
Elapsed time is 0.000185 seconds.

tic;for ii=1:N+1:M, a(ii)=inf; end,toc
Elapsed time is 0.000075 seconds.

tic; ii=1; while ii<=M, a(ii)=inf; ii=ii+N+1; end,toc
Elapsed time is 0.000053 seconds.

"Matt " <xys@whatever.com> wrote in message ex array.
> %%%%%%%%%%%%%%%%%%%%%%
> function diagtest
>
> a=rand(3000); N=length(a); M=numel(a);
>
> tic; a(1:N+1:M) = Inf; toc
> Elapsed time is 0.000185 seconds.
>
> tic;for ii=1:N+1:M, a(ii)=inf; end,toc
> Elapsed time is 0.000075 seconds.
>
> tic; ii=1; while ii<=M, a(ii)=inf; ii=ii+N+1; end,toc
> Elapsed time is 0.000053 seconds.

This test is not valid. That is why I made a copy of the matrix a in an earlier post (probably over-cautious). Apparently Matlab is smart enough when it does In-Place acceleration to see that the elements addressed are already inf (or something??). Try this:

function diagtest

a=rand(3000); N=length(a); M=numel(a);

tic; a(1:N+1:M) = Inf; toc

a=rand(3000);

tic;for ii=1:N+1:M, a(ii)=inf; end,toc

a=rand(3000);

tic; ii=1; while ii<=M, a(ii)=inf; ii=ii+N+1; end,toc




>> diagtest
Elapsed time is 0.000226 seconds.
Elapsed time is 0.000153 seconds.
Elapsed time is 0.000150 seconds.

"Matt Fig" <spamanon@yahoo.com> wrote in message <goq22p$d6c$1@fred.mathworks.com>...

> This test is not valid. That is why I made a copy of the matrix a in an earlier post (probably over-cautious). Apparently Matlab is smart enough when it does In-Place acceleration to see that the elements addressed are already inf (or something??). Try this:
__________________

Hmmm, OK.

Well, on a different matter, I'm seeing very different behavior for sparse matrices.
I'm at a loss to explain what's going on with the while-loop.

function diagtest

N=3000;
f=@(N) sprand(N,N,0.01);
%f=@(N) rand(N);

a=f(N);
M=numel(a);

tic; a(1:N+1:M) = 0; toc
%Elapsed time is 0.002220 seconds.

a(a~=0)=rand(size(a(a~=0)));

tic; for ii=1:N, a(ii,ii)=0; end,toc
%Elapsed time is 0.002273 seconds.

a(a~=0)=rand(size(a(a~=0)));


tic; ii=1; while ii<=N, a(ii,ii)=0; ii=ii+1; end,toc
%Elapsed time is 0.004197 seconds.

I never deal with sparse matrices. Maybe Tim Davis can explain what is going on. Judging by the similarities in our timings for all the other posts (under the same conditions of course) it seems our machines are similar. But here we disagree. These are my timings on your function (fairly consistently across several runs):

Elapsed time is 0.003166 seconds.
Elapsed time is 0.001761 seconds.
Elapsed time is 0.003691 seconds.

To replace diagonal of
   
    A= [1 2 3;
          4 5 6;
          7 8 9]

with inf, we can use

A(diag(true(3,1))) = inf

"Arash A Amini" wrote in message <j70c50$ls9$1@newscl01ah.mathworks.com>...
> To replace diagonal of
>
> A= [1 2 3;
> 4 5 6;
> 7 8 9]
>
> with inf, we can use
>
> A(diag(true(3,1))) = inf
===============

That will be unnecessarily expensive for larger A.

diag(true(N,1)) will generate a very large matrix, most of which will be unused.

"Matt J" wrote in message <j71hj4$efo$1@newscl01ah.mathworks.com>...

>
> That will be unnecessarily expensive for larger A.
>
> diag(true(N,1)) will generate a very large matrix, most of which will be unused.

True, so a close alternative should be used instead:

logical(speye(N))

Bruno

"Bruno Luong" <b.luong@fogale.findmycountry> wrote in message <j720h7$aj5$1@newscl01ah.mathworks.com>...
> "Matt J" wrote in message <j71hj4$efo$1@newscl01ah.mathworks.com>...
>
> >
> > That will be unnecessarily expensive for larger A.
> >
> > diag(true(N,1)) will generate a very large matrix, most of which will be unused.
>
> True, so a close alternative should be used instead:
>
> logical(speye(N))
>
> Bruno


Still seems like unnecessary overhead. Linear indexing would have to be the simplest and most direct

A(1:N+1:N^2)=.....