Compare two multiples which could overflow

[agendum97]

I have three values:

UInt64 a, b, c, d;

I need to determine if: (a*b) <=(c*d). Naturally, just doing the
multiplication and comparing would tell me the result. However, the
values could potentially overflow. For example, if:

a = UInt64.MaxValue;
b = 2;
c = UInt64.MaxValue / 3;
d = 4;

I still need to be able to conclude ab cd. I can determine if
wrapping will occur on multiplication by doing:

if (a <= 1 || b <= 1)
{
return false;
}

return ((UInt64.MaxValue / a) < b);

But it still doesn't help me determine ab cd. Anybody have any
thoughts on how I could determine this? Thanks.

[Peter Duniho]

On Fri, 03 Oct 2008 16:35:06 -0700, <ag*******@gmail.comwrote:

I have three values:

UInt64 a, b, c, d;

I need to determine if: (a*b) <=(c*d). [...]

I'm not aware of anything built into .NET that would do this directly
while handling overflow. If you were dealing with signed 64-bit values,
you could have used the Math.DivRem(), comparing a/c against d/b
(including remainder). But with unsigned 64-bit values, I think you might
as well just do the math the long way.

Pete

[Rudy Velthuis]

ag*******@gmail.com wrote:

I have three values:

UInt64 a, b, c, d;

I need to determine if: (a*b) <=(c*d). Naturally, just doing the
multiplication and comparing would tell me the result. However, the
values could potentially overflow. For example, if:

a = UInt64.MaxValue;
b = 2;
c = UInt64.MaxValue / 3;
d = 4;

I still need to be able to conclude ab cd. I can determine if
wrapping will occur on multiplication by doing:

if (a <= 1 || b <= 1)
{
return false;
}

return ((UInt64.MaxValue / a) < b);

But it still doesn't help me determine ab cd. Anybody have any
thoughts on how I could determine this? Thanks.

Well, you want to split them up a bit, into 32 bit chunks, and check
these. First some theory. Say, n = 2^32, then

a = a0 + a1*n
b = b0 + a2*n
c = c0 + c1*n
d = d0 + d1*n

(a0 + a1*n)*(b0 + b1*n) - (c0 + c1*n)*(d0 + d1*n) 0

a0*b0 + a0*b1*n + a1*b0*n + a1*b1*n*n >
c0*d0 + c0*d1*n + c1*d0*n + c1*d1*n^2 0

a0*b0 + (a0*b1 + a1*b0)*n + a1*b1*n*n >
c0*d0 + (c0*d1 + C1*d0)*n + c1*d1*n*n

Now you do how you always compare numbers, like 321 175. You first
compare the highest digits. If they are not equal, you have your
answer. If they are equal, you compare the second digits, and so on,
and so forth.

Here you do the same. if a1*b1 c1*d1, then the result is true. If
a1*b1 < c1*d1, then the result is false. And if both are equal, you
turn to comparing the next "digit", i.e. you compare
a0*b1 + a1*b0 with c0*d1 + c1*d0, etc.etc.
UInt64 a0 = a & 0xFFFFFFFF;
UInt64 a1 = a >32;
UInt64 b0 = b & 0xFFFFFFFF;
UInt64 b1 = b >32;
UInt64 c0 = c & 0xFFFFFFFF;
UInt64 c1 = c >32;
UInt64 d0 = d & 0xFFFFFFFF;
UInt64 d1 = d >32;

UInt64 x00 = a0*b0;
UInt64 x01 = c0*d0;
// this one may overflow, so split them up even more:
int x10 = (a0*b1) & 1 + (a1*b0) & 1; // lower bit
int x11 = (c0*d1) & 1 + (c1*d0) & 1; // lower bit
UInt64 x20 = a0*b1/2 + a1*b0/2; // top 63 bits
UInt64 x21 = c0*d1/2 + c1*d0/2; // top 63 bits
if (x10 = 2)
{
x20++;
x10 = 0;
}
if (x11 = 2)
{
x21++;
x11 = 0;
}
Byte x20 =
UInt64 x30 = a1*b1;
UInt64 x31 = c1*d1;

Now if the higher parts are compared, the lower bits don't matter,
except if the higher parts are equal. So you can do:

if (x30 x31)
return true;
else if (x30 < x31)
return false;
else if (x20 x21)
return true;
else if (x20 < x21)
return false;
else if (x10 x11)
return true;
else if (x10 < x11)
return false;
else if (x00 x01)
return true;
else
return false;

FWIW, this will not overflow. Each part, x00 - x31, is in the UInt64
range.

--
Rudy Velthuis http://rvelthuis.de

"The whole problem with the world is that fools and fanatics are
always so certain of themselves, but wiser people so full of
doubts." -- Bertrand Russell

[Steve Thackery]

Very helpful and informative. Thanks, Rudy.

SteveT

[=?Utf-8?B?QWdlbmR1bQ==?=]

Thanks Rudy. However, I am having difficulty understanding one part of.
Particularly:

int x10 = (a0*b1) & 1 + (a1*b0) & 1; // lower bit
int x11 = (c0*d1) & 1 + (c1*d0) & 1; // lower bit
UInt64 x20 = a0*b1/2 + a1*b0/2; // top 63 bits
UInt64 x21 = c0*d1/2 + c1*d0/2; // top 63 bits

Shouldn't a1, b1, c1, and d1 here all be multiplied by n? After all as you
have your theory above:

a = a0 + a1*n
b = b0 + a2*n
c = c0 + c1*n
d = d0 + d1*n

I guess I don't understand how the n was able to be taken out.

(a0*b1*n) & 1 + (a1*n*b0) & 1;

Because in your code what you are actually doing is:

(a0*(b >32)) & 1 + ((a >32)*b0) & 1;

And I don't see how that is helpful since doing those shifts creates some
arbitrary numbers. I mean, if you have two numbers, a and b, and you split
the numbers in binary half (aHigh, aLow, bHigh, and bLow)... you are saying:

aLow*bHigh + aHigh*bLow

And I am not sure what that is supposed to give us.

That said, in my preliminary tests it seems your code works. But I would
like to understand the logic behind it... could you expand on it a little?

Thanks


"Rudy Velthuis" wrote:

ag*******@gmail.com wrote:

I have three values:

UInt64 a, b, c, d;

I need to determine if: (a*b) <=(c*d). Naturally, just doing the
multiplication and comparing would tell me the result. However, the
values could potentially overflow. For example, if:

a = UInt64.MaxValue;
b = 2;
c = UInt64.MaxValue / 3;
d = 4;

I still need to be able to conclude ab cd. I can determine if
wrapping will occur on multiplication by doing:

if (a <= 1 || b <= 1)
{
return false;
}

return ((UInt64.MaxValue / a) < b);

But it still doesn't help me determine ab cd. Anybody have any
thoughts on how I could determine this? Thanks.

Well, you want to split them up a bit, into 32 bit chunks, and check
these. First some theory. Say, n = 2^32, then

a = a0 + a1*n
b = b0 + a2*n
c = c0 + c1*n
d = d0 + d1*n

(a0 + a1*n)*(b0 + b1*n) - (c0 + c1*n)*(d0 + d1*n) 0

a0*b0 + a0*b1*n + a1*b0*n + a1*b1*n*n >
c0*d0 + c0*d1*n + c1*d0*n + c1*d1*n^2 0

a0*b0 + (a0*b1 + a1*b0)*n + a1*b1*n*n >
c0*d0 + (c0*d1 + C1*d0)*n + c1*d1*n*n

Now you do how you always compare numbers, like 321 175. You first
compare the highest digits. If they are not equal, you have your
answer. If they are equal, you compare the second digits, and so on,
and so forth.

Here you do the same. if a1*b1 c1*d1, then the result is true. If
a1*b1 < c1*d1, then the result is false. And if both are equal, you
turn to comparing the next "digit", i.e. you compare
a0*b1 + a1*b0 with c0*d1 + c1*d0, etc.etc.
UInt64 a0 = a & 0xFFFFFFFF;
UInt64 a1 = a >32;
UInt64 b0 = b & 0xFFFFFFFF;
UInt64 b1 = b >32;
UInt64 c0 = c & 0xFFFFFFFF;
UInt64 c1 = c >32;
UInt64 d0 = d & 0xFFFFFFFF;
UInt64 d1 = d >32;

UInt64 x00 = a0*b0;
UInt64 x01 = c0*d0;
// this one may overflow, so split them up even more:
int x10 = (a0*b1) & 1 + (a1*b0) & 1; // lower bit
int x11 = (c0*d1) & 1 + (c1*d0) & 1; // lower bit
UInt64 x20 = a0*b1/2 + a1*b0/2; // top 63 bits
UInt64 x21 = c0*d1/2 + c1*d0/2; // top 63 bits
if (x10 = 2)
{
x20++;
x10 = 0;
}
if (x11 = 2)
{
x21++;
x11 = 0;
}
Byte x20 =
UInt64 x30 = a1*b1;
UInt64 x31 = c1*d1;

Now if the higher parts are compared, the lower bits don't matter,
except if the higher parts are equal. So you can do:

if (x30 x31)
return true;
else if (x30 < x31)
return false;
else if (x20 x21)
return true;
else if (x20 < x21)
return false;
else if (x10 x11)
return true;
else if (x10 < x11)
return false;
else if (x00 x01)
return true;
else
return false;

FWIW, this will not overflow. Each part, x00 - x31, is in the UInt64
range.

--
Rudy Velthuis http://rvelthuis.de

"The whole problem with the world is that fools and fanatics are
always so certain of themselves, but wiser people so full of
doubts." -- Bertrand Russell

[Rudy Velthuis]

Agendum wrote:

Thanks Rudy. However, I am having difficulty understanding one part
of. Particularly:

int x10 = (a0*b1) & 1 + (a1*b0) & 1; // lower bit
int x11 = (c0*d1) & 1 + (c1*d0) & 1; // lower bit
UInt64 x20 = a0*b1/2 + a1*b0/2; // top 63 bits
UInt64 x21 = c0*d1/2 + c1*d0/2; // top 63 bits

Shouldn't a1, b1, c1, and d1 here all be multiplied by n? After all
as you have your theory above:

a = a0 + a1*n
b = b0 + a2*n
c = c0 + c1*n
d = d0 + d1*n

I guess I don't understand how the n was able to be taken out.

Very simple. They are like numbers in base n instead of base 10. Say
the first number is u + v*n + w*n*n and the other one is x + y*n +
z*n*n. First we must compare w*n*n with z*n*n. But, if w*n*n z*n*n,
then w z, so we only have to compare w and z, and we can skip the n*n
on both sides (since n*n is positive). If they are equal, we must
compare v and y, and if they are equal too, we compare u and x. Of
course this only works because all factors (u, v, w, x, y, z) < n(*)

I had to do the little thing with the one bit since a0*b1 can't
overflow, and a1*b0 can't overflow either, but their sum can. So I took
out the low bits, and divided the rest by 2. IOW,
a0*b1/2 + a1*b0/2 can't overflow anymore, and neither can
a0*b1/2 + a1*b0/2 + 1.

(*) all numbers a0, a1, b0, b1 etc. are max. 0xFFFF, so their product
is a maximum of 0xFFFE0001, IOW always smaller than 2^32.
--
Rudy Velthuis http://rvelthuis.de

"Democracy is where you can say what you think even if you don't
think."

[Rudy Velthuis]

Agendum wrote:

Thanks Rudy. However, I am having difficulty understanding one part
of. Particularly:

int x10 = (a0*b1) & 1 + (a1*b0) & 1; // lower bit
int x11 = (c0*d1) & 1 + (c1*d0) & 1; // lower bit
UInt64 x20 = a0*b1/2 + a1*b0/2; // top 63 bits
UInt64 x21 = c0*d1/2 + c1*d0/2; // top 63 bits

Shouldn't a1, b1, c1, and d1 here all be multiplied by n? After all
as you have your theory above:

a = a0 + a1*n
b = b0 + a2*n
c = c0 + c1*n
d = d0 + d1*n

I guess I don't understand how the n was able to be taken out.

(a0*b1*n) & 1 + (a1*n*b0) & 1;

Because in your code what you are actually doing is:

(a0*(b >32)) & 1 + ((a >32)*b0) & 1;

And I don't see how that is helpful since doing those shifts creates
some arbitrary numbers.

No, they don't. A 64 bit number is generally stored in 2 32-bit dwords.
shifting a 64 bit number to the right by 32 bits leaves you with the
top dword. Comparing 64 bit numbers can first be done on the top
dwords, and then, if inconclusive (if both are the same), on the lower
ones. I use the same principle, but for a few dwords more, since you
have a maximum of 128 bits.

Say, you have the number hex 0x123456789ABCDEF0. Shifting that by 32 to
the right produces 0x12345678.

I must do the bit twiddling because the sum of a0*b1 and a1*b0 can
overflow, but a0*b1/2 + a1*b0/2 can't. This leaves however the lowest
bit unconsidered, so I take them out first, and compare them.

Addition and comparison can be done in every level: byte level, dword
level, bit level, as long as you are sure to do the addition correctly
and ad any carry bits to the next higher part.

Trust me, this is safe and sound, although I devised it on the spot and
actually tried it out in a Delphi console program first. <g>

Delphi has range checks and overflow checks built in, and can be told
to produce an error if they occur. The results are also safe. This is
how low level electronics add and compare: bit by bit. Computers add or
compare byte by byte (8 bit computers), word by word (16 bit computers)
or dword by dword (32 bit computers). How you split them up is not
really important, as long as you take care of the carry bit(s).
--
Rudy Velthuis http://rvelthuis.de

"If everything seems under control, you're just not going fast
enough." -- Mario Andretti

[Peter Morris]

Hi Rudy

Nice to see you around. I see you are still using amusing random quotes in
your signature too :-)

--
Pete
====
http://mrpmorris.blogspot.com
http://www.capableobjects.com

[Rudy Velthuis]

Peter Morris wrote:

Hi Rudy

Nice to see you around. I see you are still using amusing random
quotes in your signature too :-)

I'm still using the same newsreader, XanaNews (the best there is).

I have been reading from and writing to this group before, though.
--
Rudy Velthuis http://rvelthuis.de

"Cole's Law: Thinly sliced cabbage."

[agendum97]

On Oct 4, 6:53*am, "Rudy Velthuis" <[e-mail]:

Very simple. They are like numbers in base n instead of base 10. Say
the first number is u + v*n + w*n*n and the other one is x + y*n +
z*n*n. First we must compare w*n*n with z*n*n. But, if w*n*n z*n*n,
then w z, so we only have to compare w and z, and we can skip the n*n
on both sides (since n*n is positive). If they are equal, we must
compare v and y, and if they are equal too, we compare u and x. Of
course this only works because all factors (u, v, w, x, y, z) < n(*)

Hey Rudy,

Thanks for explaining. I studied this quite a bit, and I am pretty
sure I understand the approach you are taking here. Basically, how I
evaluated it is:

a*d <=c*b
is
(aL+aH*n)*(dL+dH*n) <=(cL+cH*n)*(bL+bH*n)
is
(aL*dL)+(aL*dH*n)+(aH*n*dL)+(aH*n*dH*n) <=(cL*bL)+(cL*bH*n)+(cH*n*bL)
+(cH*n*bH*n)

Break each side into three pieces:

1) High product: aH*dH <=cH*bH
2) Medium product: aL*dH+aH*dL <=cL*bH+cH*bL
3) Low product: aL*dL <=cL*bL

On the medium product, n was elimated simply by dividing each side by
n:

((aL*dH*n)+(aH*n*dL))/n <=((cL*bH*n)+(cH*n*bL))/n

That is how we get: aL*dH+aH*dL <=cL*dH+cH*dL

I think the only question I have left is regarding your bitshifting.
I understand why you are bitshifting (in the worst case scenario, the
sum of each medium product could overflow). However, I don't
understand why are are adding the bits and only checking for 2. For
example, you have:

int x10 = (a0*b1) & 1 + (a1*b0) & 1; // lower bit
int x11 = (c0*d1) & 1 + (c1*d0) & 1; // lower bit
UInt64 x20 = a0*b1/2 + a1*b0/2; // top 63 bits
UInt64 x21 = c0*d1/2 + c1*d0/2; // top 63 bits
if (x10 = 2)
{
x20++;
x10 = 0;
}
if (x11 = 2)
{
x21++;
x11 = 0;
}

I think the reason you are only checking for two is because it takes 2
half-bits to get one. However, doesn't this mean it is potentially
possible the medium products on one side could be 6 (0110) and 7
(0111), and then the other side 6 (0110) and 6 (0110) and the
comparison test will show they are equal when in fact 7 6... but the
bit for the 7 was lost?

Instead of checking for 2, why not add the bits onto the sum of the
divided numbers regardless of the number of bits?

int x10 = (a0*b1) & 1 + (a1*b0) & 1; // lower bit
int x11 = (c0*d1) & 1 + (c1*d0) & 1; // lower bit
UInt64 x20 = (a0*b1/2 + a1*b0/2) + x10; // top 63 bits
UInt64 x21 = (c0*d1/2 + c1*d0/2) + x11; // top 63 bits

And then skip the x10 and x11 tests entirely?

Thanks

[Rudy Velthuis]

ag*******@gmail.com wrote:

On Oct 4, 6:53*am, "Rudy Velthuis" <[e-mail]:

Very simple. They are like numbers in base n instead of base 10. Say
the first number is u + v*n + w*n*n and the other one is x + y*n +
z*n*n. First we must compare w*n*n with z*n*n. But, if w*n*n >
z*n*n, then w z, so we only have to compare w and z, and we can
skip the n*n on both sides (since n*n is positive). If they are
equal, we must compare v and y, and if they are equal too, we
compare u and x. Of course this only works because all factors (u,
v, w, x, y, z) < n(*)

Hey Rudy,

Thanks for explaining. I studied this quite a bit, and I am pretty
sure I understand the approach you are taking here. Basically, how I
evaluated it is:

a*d <=c*b
is
(aL+aH*n)*(dL+dH*n) <=(cL+cH*n)*(bL+bH*n)
is
(aL*dL)+(aL*dH*n)+(aH*n*dL)+(aH*n*dH*n) <=>
(cL*bL)+(cL*bH*n)+(cH*n*bL) +(cH*n*bH*n)

Break each side into three pieces:

1) High product: aH*dH <=cH*bH
2) Medium product: aL*dH+aH*dL <=cL*bH+cH*bL
3) Low product: aL*dL <=cL*bL

On the medium product, n was elimated simply by dividing each side by
n:

((aL*dH*n)+(aH*n*dL))/n <=((cL*bH*n)+(cH*n*bL))/n

That is how we get: aL*dH+aH*dL <=cL*dH+cH*dL

I think the only question I have left is regarding your bitshifting.
I understand why you are bitshifting (in the worst case scenario, the
sum of each medium product could overflow). However, I don't
understand why are are adding the bits and only checking for 2. For
example, you have:

Adding two bits can have 3 possible results: 0, 1, and 2 (binary: 0 1
and 10). Only 2 will not fit in one bit, and the overflow must be added
to rest of the Medium Product.
--
Rudy Velthuis http://rvelthuis.de

"The average person thinks he isn't." -- Father Larry Lorenzoni