ldiv: Make rem < 0 if a < 0.
[mes.git] / lib / mes / div.c
1 /* -*-comment-start: "//";comment-end:""-*-
2  * GNU Mes --- Maxwell Equations of Software
3  * Copyright © 2016,2017,2018,2019 Jan (janneke) Nieuwenhuizen <janneke@gnu.org>
4  * Copyright © 2019 Danny Milosavljevic <dannym@scratchpost.org>
5  *
6  * This file is part of GNU Mes.
7  *
8  * GNU Mes is free software; you can redistribute it and/or modify it
9  * under the terms of the GNU General Public License as published by
10  * the Free Software Foundation; either version 3 of the License, or (at
11  * your option) any later version.
12  *
13  * GNU Mes is distributed in the hope that it will be useful, but
14  * WITHOUT ANY WARRANTY; without even the implied warranty of
15  * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
16  * GNU General Public License for more details.
17  *
18  * You should have received a copy of the GNU General Public License
19  * along with GNU Mes.  If not, see <http://www.gnu.org/licenses/>.
20  */
21
22 #include <mes/lib.h>
23 #include <assert.h>
24 #include <stdint.h>
25 #include <limits.h>
26
27 /*#define LONG_MIN (-(~0UL>>1)-1)*/
28
29 typedef struct
30 {
31   long quot;
32   long rem;
33 } ldiv_t;
34
35 void __mesabi_div0(void)
36 {
37   eputs(" ***MES C LIB*** divide by zero numerator=");
38   eputs("\n");
39   assert(0);
40 }
41
42 /* Compare gcc: __udivdi3 */
43 unsigned long __mesabi_uldiv(unsigned long a, unsigned long b, unsigned long* remainder)
44 {
45   unsigned long tmp;
46   if (!remainder)
47     remainder = &tmp;
48   *remainder = 0;
49   if (b == 1)
50     return a;
51   else if (b == 0)
52     __mesabi_div0();
53   else
54     {
55       unsigned long x;
56       for (x = 0; a >= b; a -= b)
57         ++x;
58       *remainder = a;
59       return x;
60     }
61 }
62
63 /* Note: Rounds towards zero.
64    Maintainer: Be careful to satisfy quot * b + rem == a.
65                That means that rem can be negative. */
66 ldiv_t ldiv(long a, long b)
67 {
68   ldiv_t result;
69   int negate_result = (a < 0) ^ (b < 0);
70   assert(b != LONG_MIN);
71   if (a != LONG_MIN)
72     {
73       int negative_a = (a < 0);
74       if (negative_a)
75         a = -a;
76       if (b < 0)
77         b = -b;
78       result.quot = __mesabi_uldiv(a, b, &result.rem);
79       if (negate_result)
80         result.quot = -result.quot;
81       if (negative_a)
82         result.rem = -result.rem;
83       return result;
84     }
85   else
86     {
87       result.rem = 0;
88       if (b < 0)
89         b = -b;
90       if (b == 1)
91         {
92           result.quot = a;
93           /* Since result.quot is already negative, don't negate it again. */
94           negate_result = !negate_result;
95         }
96       else if (b == 0)
97         __mesabi_div0();
98       else
99         {
100           long x;
101           for (x = 0; a <= -b; a += b)
102             ++x;
103           result.rem = a; /* negative */
104           result.quot = x;
105         }
106       if (negate_result)
107         result.quot = -result.quot;
108       return result;
109     }
110 }