mersenneforum.org  

Go Back   mersenneforum.org > Extra Stuff > Programming

Reply
 
Thread Tools
Old 2017-06-08, 19:18   #1
ATH
Einyen
 
ATH's Avatar
 
Dec 2003
Denmark

294410 Posts
Default mulmod failing

I thought I used this well known mulmod code in quite a few programs, so I would have noticed this before: It seems to fail when a and b are large and c is small in "a*b (mod c)". By failing I mean the program is stuck in this function or it takes a very long time compared to what it should.

For example a=245673636173; b=245673636173; and c=2047; and it also "fail" when a!=b.


I thought I did a lot of tests previously with many different values a,b,c < 2^64;


Code:
uint64_t mulmod(uint64_t a, uint64_t b, uint64_t c) {
  uint64_t d; /* to hold the result of a*b mod c */
  /* calculates a*b mod c, stores result in d */
  asm ("mov %1, %%rax;"        /* put a into rax */
       "mul %2;"               /* mul a*b -> rdx:rax */
       "div %3;"               /* (a*b)/c -> quot in rax remainder in rdx */
       "mov %%rdx, %0;"        /* store result in d */
       :"=r"(d)                /* output */
       :"r"(a), "r"(b), "r"(c) /* input */
       :"%rax", "%rdx"         /* clobbered registers */
      );
  return d;
}
ATH is offline   Reply With Quote
Old 2017-06-08, 20:03   #2
danaj
 
"Dana Jacobsen"
Feb 2011
Bangkok, TH

2·3·151 Posts
Default

In my code I have this comment:

Code:
  /* GCC on a 64-bit Intel x86, help from WraithX and Wojciech Izykowski */
  /* Beware: if (a*b)/c > 2^64, there will be an FP exception */
  static inline uint64_t _mulmod(uint64_t a, uint64_t b, uint64_t n) {
    uint64_t d, dummy;                    /* d will get a*b mod c */
    asm ("mulq %3\n\t"              /* mul a*b -> rdx:rax */
         "divq %4\n\t"              /* (a*b)/c -> quot in rax remainder in rdx */
         :"=a"(dummy), "=&d"(d)     /* output */
         :"a"(a), "rm"(b), "rm"(n)  /* input */
         :"cc"                      /* mulq and divq can set conditions */
        );
    return d;
  }
It seems the exact behavior varies -- my Mac hangs, while some other machines I have exit with an FP exception.
danaj is offline   Reply With Quote
Old 2017-06-08, 20:07   #3
rogue
 
rogue's Avatar
 
"Mark"
Apr 2003
Between here and the

134428 Posts
Default

It is safest to compute a%c and b%c first. Yes, it is slower, but in the code I work with I typically need the value of a%c and b%c.
rogue is offline   Reply With Quote
Old 2017-06-08, 20:14   #4
ATH
Einyen
 
ATH's Avatar
 
Dec 2003
Denmark

27×23 Posts
Default

Thank you. You are both correct. I found my example was working when I lowered a down below the limit a*a/2047 < 2^64.

What got me confused was, that I was sure I had tested mulmod with random values up to 2^64 for all 3 variables. But now I found that old test code, and I actually did a%c and b%c before mulmod back then, and I completely forgot

Last fiddled with by ATH on 2017-06-08 at 20:15
ATH is offline   Reply With Quote
Old 2017-06-08, 22:22   #5
retina
Undefined
 
retina's Avatar
 
"The unspeakable one"
Jun 2006
My evil lair

2·32·11·29 Posts
Default

With the values you gave your code should crash with a divide-by-zero* exception. Unless the calling code has installed an exception handler or something.

Any time the high order word is larger or equal than the divisor the result from div can't fit into a single register.

245,673,636,173 x 245,673,636,173 = 60,355,535,510,463,574,085,929

60,355,535,510,463,574,085,929 / 2,047 = 29,484,873,234,227,442,152 & 785 remainder

Note that 29,484,873,234,227,442,152 cannot fit into a single 64-bit register.

This can also be detected by observing that the high order word of 60,355,535,510,463,574,085,929 is 3,271. So even before the divide is done we can know that it will fail because 3,271 >= 2,047.

* Even though the divisor is not actually zero you will still get the divide-by-zero exception whenever the result is too large.

Last fiddled with by retina on 2017-06-08 at 22:23
retina is offline   Reply With Quote
Reply

Thread Tools


Similar Threads
Thread Thread Starter Forum Replies Last Post
MISFIT uploads failing to PrimeNet Chuck PrimeNet 21 2014-02-06 15:56
Worker 1 Keeps failing zenzu88 Software 2 2012-04-10 15:16
Tried out aliqueit.exe: ggnfs failing Greebley Aliquot Sequences 35 2010-02-13 15:23
torture test failing in 1 minute ghatothkach Hardware 5 2005-01-14 09:50
Failing Torture Test.. jugbugs Hardware 12 2004-03-25 02:37

All times are UTC. The time now is 10:48.

Wed Sep 30 10:48:43 UTC 2020 up 20 days, 7:59, 0 users, load averages: 1.96, 1.60, 1.42

Powered by vBulletin® Version 3.8.11
Copyright ©2000 - 2020, Jelsoft Enterprises Ltd.

This forum has received and complied with 0 (zero) government requests for information.

Permission is granted to copy, distribute and/or modify this document under the terms of the GNU Free Documentation License, Version 1.2 or any later version published by the Free Software Foundation.
A copy of the license is included in the FAQ.