bugGNU Octave - Bugs: bug #66882, conv2 is slower for row vectors...

 
 

bug #66882: conv2 is slower for row vectors than column vectors

Submitter:  None
Submitted:  Sat 08 Mar 2025 12:20:38 AM UTC
   
 
Category:  Octave Function Severity:  3 - Normal
Priority:  5 - Normal Item Group:  Performance
Status:  In Progress Assigned to:  None
Originator Name:  Originator Email:  -email is unavailable-
Open/Closed:  * Open Release:  * 9.2.0
Release:  Operating System:  * GNU/Linux
Fixed Release:  None Planned Release:  11.1.0 (current default)
* Mandatory Fields

Add a New Comment Rich Markup
   

Jump to the original submission

Tue 11 Mar 2025 11:00:27 AM UTC, comment #9: 

Cross-linking some related topics:


Arun Giridhar <arungiridhar>
Group Member
Tue 11 Mar 2025 04:16:39 AM UTC, comment #8: 

Adding jwe to the CC list.

This is a great analysis.  I would want to put any fix deep enough such that not just users of the Octave interpreter, but users of liboctave would also get the benefit of any fix.

The only two locations that qualify are oct-convn.cc or the Fortran code.  I have a slight leaning towards modifying the C++ since Fortran seems so old and foreign that I tend to leave it alone.  Still, It would be worth understanding to what extent we could avoid creating a copy of the matrix that is going to be transposed.

Rik <rik5>
Group administrator
Mon 10 Mar 2025 10:57:16 PM UTC, comment #7: 

I made a simpler and faster test for myself:

r = ones (1, 5e4);
c = r';

tic;  x1 = conv  (r, r);  time_row_conv  = toc
tic;  x2 = conv  (c, c);  time_col_conv  = toc
tic;  x3 = conv2 (r, r);  time_row_conv2 = toc
tic;  x4 = conv2 (c, c);  time_col_conv2 = toc

assert (x1, x2'); assert (x1, x3);  assert (x1, x4');
assert (x2, x3'); assert (x2, x4);  assert (x3, x4');


Unpatched baseline:

time_row_conv = 0.195361852645874
time_col_conv = 0.1820180416107178
time_row_conv2 = 23.99320006370544
time_col_conv2 = 0.1834909915924072


The place where the actual calculation is passed from C++ to Fortran is in oct-convn.cc, which calls Fortran functions like dconv2.f and its friends. Those in turn ultimately call BLAS routines daxpy and friends.

There are two reasons for a speed difference:

  • One is the number of calls made to daxpy. It calls it once for each column, and passes the whole column in one go. For tall and skinny matrices, it is therefore very few calls, while for the short and wide case it is lots of calls (50K times more) and each time it is passed only a single element.


  • The second reason is a magnifier of the first: when daxpy is passed a lot of data in one call, it uses multiple cores, but when passed only scalars or small vectors, it stays in single-core mode. (This effect shows up clearly when using cputime instead of tic/toc).


Between the two effects, it ends up using 24 seconds instead of 0.18 seconds.

Here is some performance hackery: I am checking in dconv2.f whether the input is tall and skinny or short and wide, and I vectorize the longer dimension in the call to daxpy.

diff -r 1c0c32d8aadf liboctave/external/blas-xtra/dconv2.f
--- a/liboctave/external/blas-xtra/dconv2.f     Thu Mar 06 10:48:15 2025 -0500
+++ b/liboctave/external/blas-xtra/dconv2.f     Mon Mar 10 18:41:58 2025 -0400
@@ -39,13 +39,31 @@ c
       double precision c(ma+mb-1,na+nb-1)
       integer i,j,k
       external daxpy
-      do k = 1,na
-        do j = 1,nb
-          do i = 1,mb
-            call daxpy(ma,b(i,j),a(1,k),1,c(i,j+k-1),1)
+      if (ma + mb >= na + nb) then
+c
+c       Tall and skinny matrices: vectorize on ma.
+c
+        do k = 1,na
+          do j = 1,nb
+            do i = 1,mb
+              call daxpy(ma,b(i,j),a(1,k),1,c(i,j+k-1),1)
+            end do
           end do
         end do
-      end do
+      else
+c
+c       Short and wide matrices: vectorize on na.
+c       This currently fails "make check" because it's a hack for vectors,
+c       not for arrays with more dimensions.
+c
+        do k = 1,ma
+          do j = 1,nb
+            do i = 1,mb
+              call daxpy(na,b(i,j),a(k,1),1,c(i,j+k-1),1)
+            end do
+          end do
+        end do
+      end if
       end subroutine

       subroutine dconv2i(ma,na,a,mb,nb,b,c)


Patched:

time_row_conv = 0.1779379844665527
time_col_conv = 0.1760709285736084
time_row_conv2 = 0.1752481460571289
time_col_conv2 = 0.1754329204559326


NOTE: This hack fails "make check" because it currently only works for vectors but not for n-dimensional arrays. It's just a test hack for now.

If this sort of hackery is the right thing to do, then it would need to be done for the "outer" function in multiple Fortran files in the blas-xtra directory.

It would be much cleaner to write this check once in oct-convn.cc before it calls the corresponding Fortran function, but at the expense of constructing temporary transposed matrices. Time-memory tradeoff. (But can it be done without constructing the transpose?)

The bigger question is whether this sort of performance improvement needs to be done at all. I do not know the typical length of arrays being convolved together, so feel free to weigh in about whether the test above is representative of real use cases or not. Maybe having the user transpose it manually is the best option.

Arun Giridhar <arungiridhar>
Group Member
Mon 10 Mar 2025 10:37:23 AM UTC, comment #6: 

The value is 1e38 so 1e23 is 'epsilon' rounding error.

Better test without big values.

P = ones (1, 20000);
Q = P':
N = 1000;
tic
for ii = 1:N
    X = conv2 (P, P);
endfor
toc
tic
for ii = 1:N
    Y = conv2 (Q, Q);
endfor
toc


Anonymous
Sun 09 Mar 2025 11:42:08 PM UTC, comment #5: 

Okay, if polynomial exponential is what is desired then that is fine.  But convolution of a column vector and convolution of a row vector appear to be different operations and therefore there should be no expectation that they take the same amount of time.

Try this code


p = ones (1, 10);
x1 = 1;
for ii = 1:40
  x1 = conv2 (x1, p);
end

p = ones (10, 1);
x2 = 1;

for ii = 1:40
  x2 = conv2 (x2, p);
end

d = x1(:) - x2;
max (abs (d))


For me, I get


ans = 9.4447e+22


so the results are very different.

Rik <rik5>
Group administrator
Sun 09 Mar 2025 07:12:59 PM UTC, comment #4: 

OP code is polynomial exponential. Changing x to y is different.

Anonymous
Sun 09 Mar 2025 04:32:18 AM UTC, comment #3: 

Actually, the trouble seems to be in your benchmarking code.  You define the variable x to 1 at the start, but then use it as the output variable of the call to conv2.  There is a lot of locking/unlocking of memory and new/delete calls.

The only problem example was the third one.  If I re-write that to


p = ones (1, 10);
x = 1;
tic
for ii = 1:10000
  y = conv2 (x, p);
end
toc


then the elapsed time on my machine is 0.0271108 seconds.


Rik <rik5>
Group administrator
Sun 09 Mar 2025 01:27:40 AM UTC, comment #2: 

As Hendrik wrote, conv is a wrapper to conv2.  And, conv2 is a wrapper for convn.  And that C++ function is a wrapper to Fortran code.

What you are likely seeing is that Fortran stores matrices in column-major order.  That means that when incrementing a memory address the next matrix location fetched is one row down.  All modern processors fetch a block of memory into a cache line.  That means that the CPU is very likely to have all the row data for a given column.  However, when the for loop runs over the columns of a matrix the blocks of memory can be very far apart.  This means the data isn't in the cache and has to be fetched from DRAM over a slow memory bus.  The conv() function always makes column vectors out of its inputs so it is not a problem.

The question is whether convn, which is meant for N-dimensional objects, should have special code to detect vector inputs and re-orient them to column vectors.

Changing the Item Group to Performance and marking as In Progress.


Rik <rik5>
Group administrator
Sat 08 Mar 2025 03:54:26 AM UTC, comment #1: 

Looking at the implementation of conv will help:

type conv


As one can see, conv ALWAYS uses conv2 (covering in general 2D convolution whereas conv is a special use case for vectors only).

Note that conv converts vectors into column vectors before calling conv2.


Performance depends heavily on processor, cache size, compiler, optimization flags, used underlying specific library implementation (for conv2 I think blas is used) etc, so outside of the control of octave.

So any statement in the octave help about specific performance comparisons is not really useful.


Hendrik K <koerhen>
Sat 08 Mar 2025 12:20:38 AM UTC, original submission:  

It is very confusing whether to use conv or conv2. Somebody said "conv2 is faster than conv" because of some reason but it was actually much slower than conv but making the input a column vector suddenly becomes much faster than conv. This should be added to "help conv" and "help conv2".

Row vector and conv = Elapsed time is 3.202 seconds.

p = ones (1, 10);
x = 1;
tic
for ii = 1:10000
  x = conv (x, p);
end
toc


Column vector and conv = Elapsed time is 3.0865 seconds.

p = ones (10, 1);
x = 1;
tic
for ii = 1:10000
  x = conv (x, p);
end
toc


Row vector and conv2 = Elapsed time is 71.853 seconds.

p = ones (1, 10);
x = 1;
tic
for ii = 1:10000
  x = conv2 (x, p);
end
toc


Column vector and conv2 = Elapsed time is 1.0920 seconds.

p = ones (10, 1);
x = 1;
tic
for ii = 1:10000
  x = conv2 (x, p);
end
toc


Please, this should be added to "help conv" and "help conv2".
 
"For speed, use conv2 instead of conv, and also change row vector to column vector. Using conv2 is 3 times faster than conv for column vector but 25 times slower than conv for row vector."

Anonymous

 

(Note: upload size limit is set to 16384 kB, after insertion of the required escape characters.)

Attach Files:
   
   
Comment:
   

No files currently attached

 

Depends on the following items: None found

Items that depend on this one: None found

 

Carbon-Copy List
  • -email is unavailable- added by rik5
  • -email is unavailable- added by arungiridhar (Posted a comment)
  • -email is unavailable- added by rik5 (Posted a comment)
  • -email is unavailable- added by koerhen (Posted a comment)
  •  

    There are 0 votes so far. Votes easily highlight which items people would like to see resolved in priority, independently of the priority of the item set by tracker managers.

    Only group members can vote.

     

    Follow 6 latest changes.

    Date Changed by Updated Field Previous Value => Replaced by
    2025-03-11 rik5 Carbon-Copy- Added jwe
    2025-03-10 arungiridhar Summaryhelp conv and help conv2 should say to use column vector conv2 is slower for row vectors than column vectors
    2025-03-09 rik5 CategoryDocumentation Octave Function
        Item GroupDocumentation Performance
        StatusNone In Progress
        Planned ReleaseNone 11.1.0 (current default)

    Back to the top

    Powered by Savane 3.14-7003.
    Corresponding source code