Re: My ACCU 2016 keynote video available online
On Friday, 20 May 2016 at 20:04:35 UTC, Andrei Alexandrescu wrote: On 5/20/16 2:13 PM, Jens Müller wrote: No it doesn't work because you need to break in the last case. Consider the case when the last element of a is equal to an element in b. Next iteration you overrun a. I'm not that Bright :o). Me neither. So you'd need one more test, but you still save the other test and the test on one of the three branches, so 2 out of 4. -- Yes. Current version reduces it from 4 to 2. I've read that people use SSE to speed it up. Maybe I consider this later. If we want to improve on similar vectors (which is a great idea) we just reorder the cases, I'll guess. Just did it. It improves on my test data. But then on near dissimilar input I expect it to be worse. When doing these optimizations it is always dependent on the expected input which is a pity when optimizing library functions. These must work in all cases and should only be less efficient for very good reason. I'm looking forward to Fastware. Jens
Re: My ACCU 2016 keynote video available online
On Friday, 20 May 2016 at 14:14:18 UTC, Andrei Alexandrescu wrote:
On 05/19/2016 06:50 PM, Jens Müller wrote:
What if you stomped over an index in a that has as an equal
index in b
(it could be anywhere in b).
Hmmm, you're right. So that doesn't work, or at least not
efficiently (the fixup would entail a binary search in b).
How about this idea: arrange things such that a.back.idx <=
b.back.idx (i.e. swap them if not so). Then stomp b.back.idx
with size_t.max. Your kernel is:
if (a[i].idx < b[j].idx) {
if (i == imax) break; // check needed
i++;
} else if (a[i].idx > b[j].idx) {
j++; // no check needed
} else {
// no check needed
r += a[i++].val * b[j++].val;
}
The fixup needs only to account for the case a.back.idx ==
b.back.idx. In fact I like this a lot better than others
because it works real fast for identical vector (taking the
norm) or similar vectors (often the case). Works?
No it doesn't work because you need to break in the last case.
Consider the case when the last element of a is equal to an
element in b. Next iteration you overrun a.
But optimizing for similar vectors is interesting.
Jens
Re: My ACCU 2016 keynote video available online
On Thursday, 19 May 2016 at 22:04:56 UTC, Andrei Alexandrescu wrote: On 05/19/2016 05:36 PM, Jens Müller wrote: I removed the code to optimize for large gaps. Because it is only confusing. I may generate some benchmark data with larger gaps later to see whether it is worthwhile for such data. For skipping large gaps quickly, check galloping search (google for it, we also have it in phobos). -- Andrei Sure. I've already seen this. It's nice. But you have to include it in the sparse dot product (or list intersection) algorithm somehow. Then you require random access and galloping is only beneficial if the gaps are large. As a library writer this is a difficult position because this turns easily into over engineering. Optimally one just exposes the primitives and the user plugs them together. Ideally without having to many knobs per algorithm. Jens
Re: My ACCU 2016 keynote video available online
On Thursday, 19 May 2016 at 22:02:53 UTC, Andrei Alexandrescu
wrote:
On 05/19/2016 05:36 PM, Jens Müller wrote:
I'm not seeing it. Let me explain.
Consider the input a = [1] and b = [2, 3] (I only write the
indices).
The smallest back index is 1, i.e., a.back is the chosen
sentinel.
Nonono, you stamp the largest index over the smaller index. So
you overwrite a = [3] and you leave b = [2, 3] as it is.
Now you know that you're multiplying two correct sparse vectors
in which _definitely_ the last elements have equal indexes. So
the inner loop is:
if (a[i].idx < b[j].idx) {
i++; // no check needed
} else if (a[i].idx > b[j].idx) {
j++; // no check needed
} else {
// check needed
r += a[i].val * b[j].val;
if (i == imax || j == jmax) break;
++i;
++j;
}
At the end you need a fixup to make sure you account for the
last index that you overwrote (which of course you need to save
first).
Makes sense?
What if you stomped over an index in a that has as an equal index
in b (it could be anywhere in b). After the loop finishes you
restore the index in a. But how do you address the values for the
stomped over index if needed?
For instance test it on
a = [2]
b = [2,3]
Note the 2 in b could be anywhere.
I think you can check for
if (a[i].idx == sentinelIdx) break;
instead of
if (i == imax || j == jmax) break;
Jens
Re: My ACCU 2016 keynote video available online
On Thursday, 19 May 2016 at 12:04:31 UTC, Andrei Alexandrescu
wrote:
On 5/19/16 4:12 AM, Jens Müller wrote:
---
if (a.length == 0 || b.length == 0)
return 0;
const amax = a.length - 1, bmax = b.length - 1;
size_t i,j = 0;
double sum = 0;
for (;;)
{
if (a[i].index < b[j].index) {
if (i++ == amax) break;
}
else if (a[i].index > b[j].index) {
bumpJ: if (j++ == bmax) break;
}
else
{
assert(a[i].index == b[j].index);
sum += a[i].value * b[j].value;
if (i++ == amax) break;
goto bumpJ;
}
}
return sum;
---
Then if you add the sentinel you only need the bounds tests in
the third case.
I'm not seeing it. Let me explain.
Consider the input a = [1] and b = [2, 3] (I only write the
indices). The smallest back index is 1, i.e., a.back is the
chosen sentinel. Now I assume that we set b.back to a.back
restoring it after the loop. Now in the case a[i].index <
b[j].index I have to check whether a[i].index == a.back.index to
break because otherwise i is incremented (since a[i].index = 1
and b[j].index = 2, for i = 0 and j = 0 respectively). In the
last case I only check a[i].index == a.back.index, since this
implies b[j].index == a.back.index.
So in sum I have two bounds tests. But I think this is not what
you are thinking of.
This does not look right.
Here are the plots for the implementations.
https://www.scribd.com/doc/313204510/Running-Time
https://www.scribd.com/doc/313204526/Speedup
dot1 is my baseline, which is indeed worse than your baseline
(dot2). But only on gdc. I choose dot2 as the baseline for
computing the speedup. dot3 is the sentinel version.
I removed the code to optimize for large gaps. Because it is only
confusing. I may generate some benchmark data with larger gaps
later to see whether it is worthwhile for such data.
It looks much more regular now (ldc is still strange).
Jens
Re: My ACCU 2016 keynote video available online
On Thursday, 19 May 2016 at 12:04:31 UTC, Andrei Alexandrescu
wrote:
On 5/19/16 4:12 AM, Jens Müller wrote:
What test data did you use?
An instance for benchmarking is generated as follows. Given nnz
which is the sum of non-zero indices in input vector a and b.
auto lengthA = uniform!"[]"(0, nnz, gen);
auto lengthB = nnz - lengthA;
auto a = iota(0, nnz).randomSample(lengthA, gen).map!(i =>
Pair(i, 10)).array();
auto b = iota(0, nnz).randomSample(lengthB, gen).map!(i =>
Pair(i, 10)).array();
So I take a random sample of (0, ..., nnz) for each input.
Any better idea? I've seen that people generate vectors with
larger gaps.
10%-20% win on dot product is significant because for many
algorithms dot product is a kernel and dominates everything
else. For those any win goes straight to the bottom line.
Sure. Still I wasn't sure whether I got the idea from your talk.
So maybe there is/was more.
The base line (dot1 in the graphs) is the straightforward
version
---
size_t i,j = 0;
double sum = 0;
while (i < a.length && j < b.length)
{
if (a[i].index < b[j].index) i++;
else if (a[i].index > b[j].index) j++;
else
{
assert(a[i].index == b[j].index);
sum += a[i].value * b[j].value;
i++;
j++;
}
}
return sum;
---
That does redundant checking. There's a better baseline:
---
if (a.length == 0 || b.length == 0)
return 0;
const amax = a.length - 1, bmax = b.length - 1;
size_t i,j = 0;
double sum = 0;
for (;;)
{
if (a[i].index < b[j].index) {
if (i++ == amax) break;
}
else if (a[i].index > b[j].index) {
bumpJ: if (j++ == bmax) break;
}
else
{
assert(a[i].index == b[j].index);
sum += a[i].value * b[j].value;
if (i++ == amax) break;
goto bumpJ;
}
}
return sum;
---
I check that.
Then if you add the sentinel you only need the bounds tests in
the third case.
I post the sentinel code later. Probably there is something to
improve there as well.
BTW the effects vary greatly for different compilers.
For example with dmd the optimized version is slowest. The
baseline is
best. Weird. With gdc the optimized is best and gdc's code is
always
faster than dmd's code. With ldc it's really strange. Slower
than dmd. I
assume I'm doing something wrong here.
Used compiler flags
dmd v2.071.0
-wi -dw -g -O -inline -release -noboundscheck
gdc (crosstool-NG 203be35 - 20160205-2.066.1-e95a735b97) 5.2.0
-Wall -g -O3 -fomit-frame-pointer -finline-functions -frelease
-fno-bounds-check -ffast-math
ldc (0.15.2-beta2) based on DMD v2.066.1 and LLVM 3.6.1
-wi -dw -g -O3 -enable-inlining -release -boundscheck=off
Am I missing some flags?
These look reasonable.
But ldc looks so bad.
Any comments from ldc users or developers? Because I see this in
many other measurements as well. I would love to have another
compiler producing efficient like gdc.
For example what's equivalent to gdc's -ffast-math in ldc.
I uploaded my plots.
- running time
https://www.scribd.com/doc/312951947/Running-Time
- speed up https://www.scribd.com/doc/312951964/Speedup
What is dot2? Could you please redo the experiments with the
modified code as well?
dot2 is an optimization for jumping over gaps more quickly
replacing the first two if statements with while statements.
But my benchmark tests have no large gaps but interestingly it
does make things worse.
Jens
Re: My ACCU 2016 keynote video available online
On Monday, 16 May 2016 at 13:46:11 UTC, Andrei Alexandrescu wrote:
Uses D for examples, showcases Design by Introspection, and
rediscovers a fast partition routine. It was quite well
received. https://www.youtube.com/watch?v=AxnotgLql0k
Andrei
Nice presentation.
The code applying the sentinel optimization assumes mutability of
the input.
That needs to be checked for. That's fine for partition because
that is assumed
to be in-place. But for other algorithms it's not so obvious.
It's sad that the
optimization works only for non-const input. It is in conflict
with the advice
to make input const if the function doesn't change it. This makes
the
optimization less likely to be applicable. One might though relax
the const
requirement to mean "the input is identical at return of the
function to its
beginning". But that's a different story, I'll guess. Coming up
with another
implementation might also work, using chain or so. But typically
the sentinel
optimization assumes mutability.
I didn't get the idea behind sentinels for sparse dot product. I
picked the
smallest of the last elements (so you need bidirectional ranges)
and fix up as
needed. For gdc I get a speedup (baseline over new
implementation) of 1.2 in
best case and >1.0 in worst case. On average it's about 1.1 I
would say. I
expected more. How would you approach sentinels with the sparse
dot product. Can
you elaborate the idea from the video? I didn't get it.
The base line (dot1 in the graphs) is the straightforward version
---
size_t i,j = 0;
double sum = 0;
while (i < a.length && j < b.length)
{
if (a[i].index < b[j].index) i++;
else if (a[i].index > b[j].index) j++;
else
{
assert(a[i].index == b[j].index);
sum += a[i].value * b[j].value;
i++;
j++;
}
}
return sum;
---
BTW the effects vary greatly for different compilers.
For example with dmd the optimized version is slowest. The
baseline is
best. Weird. With gdc the optimized is best and gdc's code is
always
faster than dmd's code. With ldc it's really strange. Slower than
dmd. I
assume I'm doing something wrong here.
Used compiler flags
dmd v2.071.0
-wi -dw -g -O -inline -release -noboundscheck
gdc (crosstool-NG 203be35 - 20160205-2.066.1-e95a735b97) 5.2.0
-Wall -g -O3 -fomit-frame-pointer -finline-functions -frelease
-fno-bounds-check -ffast-math
ldc (0.15.2-beta2) based on DMD v2.066.1 and LLVM 3.6.1
-wi -dw -g -O3 -enable-inlining -release -boundscheck=off
Am I missing some flags?
I uploaded my plots.
- running time https://www.scribd.com/doc/312951947/Running-Time
- speed up https://www.scribd.com/doc/312951964/Speedup
*Disclaimer*
I hope most of this makes sense but take it with a grain of salt.
Jens
PS
It seems the mailinglist interface does not work. I cannot send
replies anymore via mail. I wrote Brad Roberts but no answer yet.
