source: project/release/4/numbers/trunk/numbers-c.c @ 31291

Last change on this file since 31291 was 31291, checked in by sjamaan, 4 years ago

numbers: convert negation of bignums to core naming convention. This includes an on-stack/native heap C_allocate_bignum (in CPS context)

File size: 73.2 KB
Line 
1/* numbers-c.c
2 *
3 * Copyright 2010-2014 The CHICKEN Team
4 *
5 * This contains a barely recognizable version of c/bignum.c from Scheme48 1.8:
6 * Copyright (c) 1993-2008 Richard Kelsey and Jonathan Rees
7 * Copyright 1986,1987,1988,1989,1990,1991 Massachusetts Institute of Technology
8 * Copyright 1992,1993,1994,2004 Massachusetts Institute of Technology
9 *
10 * Redistribution and use in source and binary forms, with or without
11 * modification, are permitted provided that the following conditions are
12 * met:
13 *
14 *    1. Redistributions of source code must retain the above copyright
15 *       notice, this list of conditions and the following disclaimer.
16 *
17 *    2. Redistributions in binary form must reproduce the above
18 *       copyright notice, this list of conditions and the following
19 *       disclaimer in the documentation and/or other materials provided
20 *       with the distribution.
21 *
22 *    3. The name of the author may not be used to endorse or promote
23 *       products derived from this software without specific prior
24 *       written permission.
25 *
26 * THIS SOFTWARE IS PROVIDED BY THE AUTHOR ``AS IS'' AND ANY EXPRESS OR
27 * IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE IMPLIED
28 * WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE ARE
29 * DISCLAIMED. IN NO EVENT SHALL THE AUTHOR BE LIABLE FOR ANY DIRECT,
30 * INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES
31 * (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR
32 * SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION)
33 * HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT,
34 * STRICT LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING
35 * IN ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE
36 * POSSIBILITY OF SUCH DAMAGE.
37 */
38
39#include <assert.h>
40#include <errno.h>
41
42static void *tags;
43
44#include "numbers-c.h"
45
46#define fix_to_flo(p, n, f)       C_flonum(p, C_unfix(f))
47#define big_of(v)                 ((bignum_type)C_data_pointer(C_block_item(v, 1)))
48
49static C_word
50init_tags(___scheme_value tagvec)
51{
52  tags = CHICKEN_new_gc_root();
53  CHICKEN_gc_root_set(tags, tagvec);
54  return C_SCHEME_UNDEFINED;
55}
56
57/*
58 * This is an odd one out. It doesn't accept a continuation.
59 * I've put in the (more or less) verbatim code for s48_bignum_to_double.
60 */
61static C_word
62big_to_flo(C_word **p, C_word n, C_word value)
63{
64  bignum_type bignum = big_of(value);
65  double accumulator = 0;
66  bignum_digit_type * start = (BIGNUM_START_PTR (bignum));
67  bignum_digit_type * scan = (start + (BIGNUM_LENGTH (bignum)));
68  while (start < scan)
69    accumulator = ((accumulator * BIGNUM_RADIX) + (*--scan));
70  return C_flonum(p, (BIGNUM_NEGATIVE_P(bignum) ? -accumulator : accumulator));
71}
72
73
74static C_word
75check_number(C_word x)
76{
77  C_word tagvec;
78
79  if((x & C_FIXNUM_BIT) != 0) return C_fix(FIX);
80  else if(C_immediatep(x)) return C_fix(NONE);
81  else if(C_header_bits(x) == C_FLONUM_TYPE) return C_fix(FLO);
82  else if(C_header_bits(x) == C_STRUCTURE_TYPE) {
83    tagvec = CHICKEN_gc_root_ref(tags);
84   
85    if (C_block_item(x, 0) == C_block_item(tagvec, BIG_TAG))
86      return C_fix(BIG);
87    else if (C_block_item(x, 0) == C_block_item(tagvec, COMP_TAG))
88      return C_fix(COMP);
89    else if (C_block_item(x, 0) == C_block_item(tagvec, RAT_TAG))
90      return C_fix(RAT);
91    else
92      return C_fix(NONE);
93  } else
94      return C_fix(NONE);
95}
96
97static void
98C_wrap_bignum(C_word c, C_word closure, C_word bigvec)
99{
100  C_word ab[3], *a = ab;
101  bignum_type src = (bignum_type)C_block_item(C_block_item(closure, 2), 0);
102  bignum_digit_type *scan = BIGNUM_TO_POINTER(src);
103  bignum_digit_type *end = BIGNUM_START_PTR(src) + BIGNUM_LENGTH(src);
104  bignum_digit_type *tgt = C_data_pointer(bigvec);
105  C_word tagvec = CHICKEN_gc_root_ref(tags);
106  C_word result = C_structure(&a, 2, C_block_item(tagvec, BIG_TAG), bigvec);
107  while (scan < end)
108    (*tgt++) = (*scan++);
109  BIGNUM_DEALLOCATE(src);
110  C_kontinue(C_block_item(closure, 1), result);
111}
112
113static void
114C_return_bignum(C_word k, bignum_type b)
115{
116  C_word result;
117 
118  assert(b != BIGNUM_OUT_OF_BAND);
119
120  switch(BIGNUM_LENGTH(b)) {
121  case 0:
122    BIGNUM_DEALLOCATE(b);
123    C_kontinue(k, C_fix(0));
124  /* This code "knows" that bignums have 2 "reserved" bits, like fixnums */
125  case 1:
126    result = C_fix(BIGNUM_NEGATIVE_P(b) ? -BIGNUM_REF(b, 0) : BIGNUM_REF(b, 0));
127    BIGNUM_DEALLOCATE(b);
128    C_kontinue(k, result);
129  case 2:
130    /* Edge case: most negative fixnum */
131    if (BIGNUM_NEGATIVE_P(b) && BIGNUM_REF(b, 0) == 0 && BIGNUM_REF(b, 1) == 1) {
132      BIGNUM_DEALLOCATE(b);
133      C_kontinue(k, C_fix(C_MOST_NEGATIVE_FIXNUM));
134    }
135    /* FALLTHROUGH */
136  default:
137    {
138      C_word pab[2], *pa = pab, kab[4], *ka = kab, k2;
139      /* Make result a wrapped pointer because C_closure wants scheme objects */
140      result = C_mpointer(&pa, b);
141      k2 = C_closure(&ka, 3, (C_word)C_wrap_bignum, k, result);
142      /* Here we assume bignum digits are C words.. */
143      C_allocate_vector(6, (C_word)NULL, k2,
144                        C_fix(sizeof(C_word) * (BIGNUM_LENGTH(b) + 1)),
145                        /* Byte vec, no initialization, align at 8 bytes */
146                        C_SCHEME_TRUE, C_SCHEME_FALSE, C_SCHEME_FALSE);
147    }
148  }
149}
150
151static void
152C_bignum_wrapped_return_bigobj(C_word c, C_word closure, C_word wrapped_big)
153{
154  C_word obj = C_block_item(closure, 2);
155  C_values(4, C_SCHEME_UNDEFINED, C_block_item(closure, 1), wrapped_big, obj);
156}
157
158static void
159C_return_big_fix(C_word k, bignum_type big, C_word fix)
160{
161  C_word kab[4], *ka = kab, k2;
162  k2 = C_closure(&ka, 3, (C_word)C_bignum_wrapped_return_bigobj, k, fix);
163  C_return_bignum(k2, big);
164}
165
166static void
167C_b2_wrapped(C_word c, C_word closure, C_word wrapped_b2)
168{
169  C_word k = C_block_item(closure, 1); /* Original closure */
170  C_word kab[4], *ka = kab, k2; /* Next continuation */
171  k2 = C_closure(&ka, 3, (C_word)C_bignum_wrapped_return_bigobj, k, wrapped_b2);
172  C_return_bignum(k2, (bignum_type)C_block_item(C_block_item(closure, 2), 0));
173}
174
175static void
176C_return_2_bignums(C_word k, bignum_type b1, bignum_type b2)
177{
178  C_word bab[2], *ba = bab, kab[4], *ka = kab, k2, b1_ptr;
179  /* Make b1 a wrapped pointer because C_closure wants scheme objects */
180  b1_ptr = C_mpointer(&ba, b1);
181  /* Wrap b2 first, then b1. Return them to k in the original (b1,b2) order */
182  k2 = C_closure(&ka, 3, (C_word)C_b2_wrapped, k, b1_ptr);
183  C_return_bignum(k2, b2);
184}
185
186
187static bignum_type
188bignum_allocate(bignum_length_type length, int negative_p)
189{
190  bignum_type result;
191  bignum_digit_type *digits;
192
193  digits = (bignum_digit_type *)C_malloc(sizeof(bignum_digit_type)*(length + 1));
194 
195  if(digits == NULL) {
196    fprintf(stderr, "out of memory - can not allocate number");
197    exit(EXIT_FAILURE);
198  }
199
200  result = (bignum_type)digits;
201  BIGNUM_SET_HEADER(result, length, negative_p);
202  return result;
203}
204
205static bignum_type
206bignum_allocate_zeroed(bignum_length_type length, int negative_p)
207{
208  assert((length >= 0) || (length < BIGNUM_RADIX));
209  {
210    bignum_type result;
211    bignum_digit_type *digits;
212    bignum_digit_type *scan;
213    bignum_digit_type *end;
214    digits=(bignum_digit_type *)C_malloc(sizeof(bignum_digit_type)*(length+1));
215 
216    if(digits == NULL) {
217      fprintf(stderr, "out of memory - can not allocate number");
218      exit(EXIT_FAILURE);
219    }
220
221    result = (bignum_type)digits;
222    BIGNUM_SET_HEADER(result, length, negative_p);
223    scan = BIGNUM_START_PTR(result);
224    end = scan + length;
225    while (scan < end)
226      (*scan++) = 0;
227    return (result);
228  }
229}
230
231static bignum_type
232bignum_allocate_from_fixnum(C_word fix)
233{
234  bignum_type ret;
235
236  if (fix == C_fix(0)) {
237    ret = bignum_allocate(0, 0);
238  } else if (fix == C_fix(C_MOST_NEGATIVE_FIXNUM)) {
239    ret = bignum_allocate(2, 1);
240    BIGNUM_REF(ret, 0) = 0;
241    BIGNUM_REF(ret, 1) = 1;
242  } else {
243    bignum_digit_type digit = C_unfix(fix);
244    ret = bignum_allocate(1, digit < 0);
245    BIGNUM_REF(ret, 0) = ((digit < 0) ? -digit : digit);
246  }
247  return ret;
248}
249
250static bignum_type
251bignum_digit_to_bignum(bignum_digit_type digit, int neg_p)
252{
253  bignum_type ret;
254  if (digit == 0)
255    return bignum_allocate(0, 0);
256
257  ret = bignum_allocate(1, neg_p);
258  BIGNUM_REF(ret, 0) = digit;
259  return ret;
260}
261
262static bignum_type
263shorten_bignum(bignum_type big, bignum_length_type newlength)
264{
265  bignum_digit_type *digits, *newdigits;
266  digits = BIGNUM_TO_POINTER(big);
267  newdigits = (bignum_digit_type *)C_realloc(digits, sizeof(bignum_digit_type)*(newlength + 1));
268  if (newdigits == NULL) {
269    fprintf(stderr, "out of memory - can not reallocate number");
270    exit(EXIT_FAILURE);
271  }
272  return (bignum_type)newdigits;
273}
274
275static bignum_type
276bignum_trim(bignum_type bignum)
277{
278  bignum_digit_type * start = (BIGNUM_START_PTR (bignum));
279  bignum_digit_type * end = (start + (BIGNUM_LENGTH (bignum)));
280  bignum_digit_type * scan = end;
281  while ((start <= scan) && ((*--scan) == 0))
282    ;
283  scan += 1;
284  if (scan < end)
285    {
286      bignum_length_type length = (scan - start);
287      BIGNUM_SET_HEADER
288        (bignum, length, ((length != 0) && (BIGNUM_NEGATIVE_P (bignum))));
289      BIGNUM_REDUCE_LENGTH (bignum, bignum, length);
290    }
291  return (bignum);
292}
293
294/* Copying */
295
296static bignum_type
297bignum_copy(bignum_type source)
298{
299  bignum_type target =
300    (bignum_allocate ((BIGNUM_LENGTH (source)), (BIGNUM_NEGATIVE_P (source))));
301  bignum_destructive_copy (source, target);
302  return (target);
303}
304
305static void
306bignum_destructive_copy(bignum_type source, bignum_type target)
307{
308  bignum_digit_type * scan_source = (BIGNUM_START_PTR (source));
309  bignum_digit_type * end_source =
310    (scan_source + (BIGNUM_LENGTH (source)));
311  bignum_digit_type * scan_target = (BIGNUM_START_PTR (target));
312  while (scan_source < end_source)
313    (*scan_target++) = (*scan_source++);
314  return;
315}
316
317static bignum_type
318bignum_new_sign(bignum_type bignum, int negative_p)
319{
320  bignum_type result =
321    (bignum_allocate ((BIGNUM_LENGTH (bignum)), negative_p));
322  bignum_destructive_copy (bignum, result);
323  return (result);
324}
325
326#define fix_randomize  C_randomize
327#define fix_random C_random_fixnum
328
329/*
330 * This random number generator is very simple. Probably too simple...
331 */
332static C_word
333big_randomize(C_word bignum)
334{
335  bignum_type big;
336  C_word seed = 0;
337  bignum_digit_type * scan, * end;
338
339  big = big_of(bignum);
340  scan = (BIGNUM_START_PTR (big));
341  end = scan + BIGNUM_LENGTH(big);
342  /* What a cheap way to initialize the random generator. I feel dirty! */
343  while (scan < end)
344    seed ^= *scan++;
345
346  srand(C_unfix(seed));
347  return C_SCHEME_UNDEFINED;
348}
349
350static void
351big_random(C_word c, C_word self, C_word k, C_word max)
352{
353  bignum_type bigmax;
354  bignum_type result;
355  bignum_length_type randlen, max_len, max_bits;
356  bignum_digit_type max_top_digit, d;
357  bignum_digit_type * scan, * end;
358 
359  bigmax = big_of(max);
360
361  max_len = BIGNUM_LENGTH(bigmax);
362  max_top_digit = d = BIGNUM_REF(bigmax, max_len - 1);
363 
364  max_bits = (max_len - 1) * BIGNUM_DIGIT_LENGTH;
365  while(d) {
366    max_bits++;
367    d >>= 1;
368  }
369  /* Subtract/add one because we don't want zero to be over-represented */
370  randlen = (((double)rand())/(RAND_MAX + 1.0) * (double)(max_bits - 1));
371  randlen = BIGNUM_BITS_TO_DIGITS(randlen + 1);
372 
373  result = bignum_allocate(randlen, BIGNUM_NEGATIVE_P(bigmax));
374  scan = BIGNUM_START_PTR(result);
375  end = scan + randlen - 1; /* Go to just before the end */
376  while(scan < end)
377    *scan++ = (((double)rand())/(RAND_MAX + 1.0) * (double)BIGNUM_RADIX);
378  /*
379   * Last word is special when length is max_len: It must be less than
380   * max's most significant digit, instead of BIGNUM_RADIX.
381   */
382  if (max_len == randlen)
383    *scan = (((double)rand())/(RAND_MAX + 1.0) * (double)max_top_digit);
384  else
385    *scan = (((double)rand())/(RAND_MAX + 1.0) * (double)BIGNUM_RADIX);
386
387  result = bignum_trim(result); 
388  C_return_bignum(k, result);
389}
390
391static void
392fix_plus_big(C_word c, C_word self, C_word k, C_word x, C_word y)
393{
394  bignum_type bigx, result;
395  bigx = bignum_allocate_from_fixnum(x);
396  result = bignum_add(bigx, big_of(y));
397  BIGNUM_DEALLOCATE(bigx);
398  C_return_bignum(k, result);
399}
400
401
402static void
403fix_minus_big(C_word c, C_word self, C_word k, C_word x, C_word y)
404{
405  bignum_type tmp, result;
406  tmp = bignum_allocate_from_fixnum(x);
407
408  result = bignum_subtract(tmp, big_of(y));
409  BIGNUM_DEALLOCATE(tmp);
410  C_return_bignum(k, result);
411}
412
413
414static void
415big_minus_fix(C_word c, C_word self, C_word k, C_word x, C_word y)
416{
417  bignum_type tmp, result;
418  tmp = bignum_allocate_from_fixnum(y);
419
420  result = bignum_subtract(big_of(x), tmp);
421  BIGNUM_DEALLOCATE(tmp);
422  C_return_bignum(k, result);
423}
424
425static void
426big_divrem_big(C_word c, C_word self, C_word k, C_word x, C_word y)
427{
428  bignum_type numerator = big_of(x);
429  bignum_type denominator = big_of(y);
430  int q_neg_p = ((BIGNUM_NEGATIVE_P (denominator))
431                 ? (! (BIGNUM_NEGATIVE_P (numerator)))
432                 : (BIGNUM_NEGATIVE_P (numerator)));
433
434  switch (bignum_compare_unsigned (numerator, denominator))
435    {
436    case bignum_comparison_equal:
437      C_values(4, C_SCHEME_UNDEFINED, k,
438               q_neg_p ? C_fix(-1) : C_fix(1), C_fix(0));
439    case bignum_comparison_less:
440      C_values(4, C_SCHEME_UNDEFINED, k, C_fix(0), x);
441    case bignum_comparison_greater:
442    default:                                  /* to appease gcc -Wall */
443      {
444        bignum_type quotient, remainder;
445        int r_neg_p = BIGNUM_NEGATIVE_P (numerator) ? 1 : 0;
446       
447        assert(BIGNUM_LENGTH(denominator) > 1);
448        bignum_divide_unsigned_large_denominator
449          (numerator, denominator, (&quotient), (&remainder), q_neg_p, r_neg_p);
450        C_return_2_bignums(k, quotient, remainder);
451      }
452    }
453}
454
455static void
456big_divrem_fix(C_word c, C_word self, C_word k, C_word x, C_word y)
457{
458  bignum_type bigx = big_of(x);
459  y = C_unfix(y);
460
461  if (y == 1)
462    C_values(4, C_SCHEME_UNDEFINED, k, x, C_fix(0));
463  else if (y == -1)
464    C_return_big_fix(k, bignum_new_sign(bigx, !(BIGNUM_NEGATIVE_P(bigx))), C_fix(0));
465
466  /* Too bad, we really need to do some work... */
467  {
468    int q_neg_p = (y < 0) ? !(BIGNUM_NEGATIVE_P(bigx)) : BIGNUM_NEGATIVE_P(bigx);
469    int r_neg_p = BIGNUM_NEGATIVE_P(bigx);
470    bignum_digit_type abs_y = (y < 0) ? -y : y;
471    bignum_type quotient, remainder;
472   
473    if (y == C_MOST_NEGATIVE_FIXNUM) {
474      if (!BIGNUM_NEGATIVE_P(bigx) && BIGNUM_LENGTH(bigx) == 1
475          && BIGNUM_REF(bigx, 1) == 1 && BIGNUM_REF(bigx, 0) == 0) {
476        /*
477         * Very very special case:
478         * quotient(MOST_NEGATIVE_FIXNUM, -(MOST_NEGATIVE_FIXNUM)) => -1
479         */
480        C_values(4, C_SCHEME_UNDEFINED, k, C_fix(-1), C_fix(0));
481      } else {
482        /* This is the only case we need to go allocate a bignum for */
483        bignum_type bigy =
484          bignum_allocate_from_fixnum(C_fix(C_MOST_NEGATIVE_FIXNUM));
485
486        bignum_divide_unsigned_large_denominator
487          (bigx, bigy, (&quotient), (&remainder), q_neg_p, r_neg_p);
488        BIGNUM_DEALLOCATE(bigy);
489        C_return_2_bignums(k, quotient, remainder);
490      }
491    } else if (abs_y < BIGNUM_RADIX_ROOT) {
492      bignum_divide_unsigned_small_denominator
493        (bigx, abs_y, (&quotient), (&remainder), q_neg_p, r_neg_p);
494      C_return_2_bignums(k, quotient, remainder);
495    } else {
496      bignum_divide_unsigned_medium_denominator
497        (bigx, abs_y, (&quotient), (&remainder), q_neg_p, r_neg_p);
498      C_return_2_bignums(k, quotient, remainder);
499    }
500  }
501}
502
503/*
504 * Big_gcd is a huge function and it sucks that it needs to be in C.
505 *
506 * Why does it need to be in C?  Because if you have a very big bignum
507 * that doesn't cleanly divide another big bignum, you end up calling
508 * the remainder procedure a lot in Scheme.  This produces tons of
509 * intermediate bignums, which means a lot of copies into GC'able memory
510 * need to be made (and the GC will be triggered more often).
511 * That's a major slowdown.  Doing the loop in C means the intermediate
512 * results can be cleaned up right away each loop step, and returning
513 * just one result to Scheme.
514 * Once (if?) we find a way to avoid copying bignums (instead allocating
515 * directly in GCable memory) this function can be cut out and replaced by
516 * a recursive call to gcd-0 in Scheme.
517 */
518static void
519big_gcd(C_word c, C_word self, C_word k, C_word x, C_word y)
520{
521  bignum_type bigx = big_of(x), bigy = big_of(y), bigr;
522
523  assert(BIGNUM_LENGTH(bigx) > 1);
524  assert(BIGNUM_LENGTH(bigy) > 1);
525 
526  switch(bignum_compare_unsigned (bigx, bigy)) {
527  case bignum_comparison_equal:
528    C_kontinue(k, x);
529  case bignum_comparison_less:
530    /* Swap, since remainder of bigx, bigy would be bigx, causing an extra loop */
531    bigr = bigy;
532    bigy = bigx;
533    bigx = bigr;
534
535    /* Make x and y match for the special case where gcd(x, y) = y */
536    {
537      C_word tmp = y;
538      y = x;
539      x = tmp;
540    }
541    /* FALLTHROUGH */
542  default: /* Continue below */
543    break;
544  }
545 
546  /*
547   * Be careful! Don't deallocate live objects. We could start with a copy
548   * or compare pointers with big_of(x) or y every time but that seems wasteful.
549   */
550  bignum_divide_unsigned_large_denominator
551   (bigx, bigy, ((bignum_type *) 0), (&bigr), 0, 0);
552  bigx = bigy;
553  bigy = bigr;
554  /* Original bigx is forgotten now */
555  assert(bigx != big_of(x));
556  assert(bigy != big_of(x));
557  /* Only bigx points to y */
558  assert(bigy != big_of(y));
559  assert(bigx == big_of(y));
560
561  switch (BIGNUM_LENGTH(bigy)) {
562  case 0:  /* bigy = 0 */
563    /* remainder(x, y) = 0 => y  |  x < y */
564    BIGNUM_DEALLOCATE(bigy); /* Got allocated in previous step */
565    C_kontinue(k, y);
566  case 1:
567    if (BIGNUM_REF(bigy, 0) == 1) { /* y = 1? Then 1 is the result */
568      BIGNUM_DEALLOCATE(bigy); /* Got allocated in previous step */
569      C_kontinue(k, C_fix(1));
570    } else if (BIGNUM_REF(bigy, 0) < BIGNUM_RADIX_ROOT)
571      bigr = bignum_remainder_unsigned_small_denominator
572              (bigx, BIGNUM_REF(bigy, 0), 0);
573    else
574      bignum_divide_unsigned_medium_denominator
575       (bigx, BIGNUM_REF(bigy, 0), (bignum_type *)0, (&bigr), 0, 0);
576    break;
577  default:
578    bignum_divide_unsigned_large_denominator
579     (bigx, bigy, ((bignum_type *) 0), (&bigr), 0, 0);
580  }
581  /* Swap, but don't deallocate x since it holds the original value of y */
582  bigx = bigy;
583  bigy = bigr;
584 
585  /* Original bigy is forgotten now, we can safely always deallocate bigx */
586
587  /* Assume that bignums coming from outside are never length 1 */
588  assert(bigx != big_of(y));
589 
590  while(BIGNUM_LENGTH(bigy) > 1) {
591    bignum_divide_unsigned_large_denominator
592     (bigx, bigy, ((bignum_type *) 0), (&bigr), 0, 0);
593    BIGNUM_DEALLOCATE(bigx);
594    bigx = bigy;
595    bigy = bigr;
596  }
597
598  /* Finish up with a faster loop until y = 0 (ie, length(bigy) = 0) */
599  while (BIGNUM_LENGTH(bigy) == 1) {
600    if (BIGNUM_REF(bigy, 0) == 1) {
601      BIGNUM_DEALLOCATE(bigx);
602      BIGNUM_DEALLOCATE(bigy);
603      C_kontinue(k, C_fix(1));
604      break;
605    } else if (BIGNUM_REF(bigy, 0) < BIGNUM_RADIX_ROOT) {
606      bigr = bignum_remainder_unsigned_small_denominator
607              (bigx, BIGNUM_REF(bigy, 0), 0);
608    } else {
609      bignum_divide_unsigned_medium_denominator
610       (bigx, BIGNUM_REF(bigy, 0), (bignum_type *)0, (&bigr), 0, 0);
611    }
612    BIGNUM_DEALLOCATE(bigx);
613    bigx = bigy;
614    bigy = bigr;
615  }
616  BIGNUM_DEALLOCATE(bigy);
617  C_return_bignum(k, bigx);
618}
619
620static C_word
621C_a_i_flonum_gcd(C_word **p, C_word n, C_word x, C_word y)
622{
623   double xub, yub, r;
624   
625   xub = C_flonum_magnitude(x);
626   yub = C_flonum_magnitude(y);
627
628   if (xub < 0.0) xub = -xub;
629   if (yub < 0.0) yub = -yub;
630   
631   while(yub != 0.0) {
632     r = fmod(xub, yub);
633     xub = yub;
634     yub = r;
635   }
636   C_return(C_flonum(p, xub));
637}
638
639
640static void
641big_plus_big(C_word c, C_word self, C_word k, C_word x, C_word y)
642{
643  bignum_type result;
644  result = bignum_add(big_of(x), big_of(y));
645  C_return_bignum(k, result);
646}
647
648/*
649 * This now makes the assumption it is never passed a bignum of LENGTH 0.
650 * This should always be valid in Chicken.
651 */
652static bignum_type
653bignum_add(bignum_type x, bignum_type y)
654{
655  return
656    (((BIGNUM_NEGATIVE_P (x))
657      ? ((BIGNUM_NEGATIVE_P (y))
658         ? (bignum_add_unsigned (x, y, 1))
659         : (bignum_subtract_unsigned (y, x)))
660      : ((BIGNUM_NEGATIVE_P (y))
661         ? (bignum_subtract_unsigned (x, y))
662         : (bignum_add_unsigned (x, y, 0)))));
663}
664
665static bignum_type
666bignum_add_unsigned(bignum_type x, bignum_type y, int negative_p)
667{
668  if ((BIGNUM_LENGTH (y)) > (BIGNUM_LENGTH (x)))
669    {
670      bignum_type z = x;
671      x = y;
672      y = z;
673    }
674  {
675    bignum_length_type x_length = (BIGNUM_LENGTH (x));
676    bignum_type r = (bignum_allocate ((x_length + 1), negative_p));
677    bignum_digit_type sum;
678    bignum_digit_type carry = 0;
679    bignum_digit_type * scan_x = (BIGNUM_START_PTR (x));
680    bignum_digit_type * scan_r = (BIGNUM_START_PTR (r));
681    {
682      bignum_digit_type * scan_y = (BIGNUM_START_PTR (y));
683      bignum_digit_type * end_y = (scan_y + (BIGNUM_LENGTH (y)));
684      while (scan_y < end_y)
685        {
686          sum = ((*scan_x++) + (*scan_y++) + carry);
687          if (sum < BIGNUM_RADIX)
688            {
689              (*scan_r++) = sum;
690              carry = 0;
691            }
692          else
693            {
694              (*scan_r++) = (sum - BIGNUM_RADIX);
695              carry = 1;
696            }
697        }
698    }
699    {
700      bignum_digit_type * end_x = ((BIGNUM_START_PTR (x)) + x_length);
701      if (carry != 0)
702        while (scan_x < end_x)
703          {
704            sum = ((*scan_x++) + 1);
705            if (sum < BIGNUM_RADIX)
706              {
707                (*scan_r++) = sum;
708                carry = 0;
709                break;
710              }
711            else
712              (*scan_r++) = (sum - BIGNUM_RADIX);
713          }
714      while (scan_x < end_x)
715        (*scan_r++) = (*scan_x++);
716    }
717    if (carry != 0) {
718      (*scan_r) = 1;
719    } else { /* r is one word too big (to hold a possible carry), readjust */
720      BIGNUM_SET_HEADER(r, x_length, negative_p);
721      BIGNUM_REDUCE_LENGTH(r, r, x_length);
722    }
723    return (r);
724  }
725}
726
727static void
728big_minus_big(C_word c, C_word self, C_word k, C_word x, C_word y)
729{
730  bignum_type result;
731  result = bignum_subtract(big_of(x), big_of(y));
732  C_return_bignum(k, result);
733}
734
735/*
736 * This now makes the assumption it is never passed a bignum of LENGTH 0.
737 * This should always be valid in Chicken.
738 */
739static bignum_type
740bignum_subtract(bignum_type x, bignum_type y)
741{
742  return
743    (((BIGNUM_NEGATIVE_P (x))
744       ? ((BIGNUM_NEGATIVE_P (y))
745          ? (bignum_subtract_unsigned (y, x))
746          : (bignum_add_unsigned (x, y, 1)))
747       : ((BIGNUM_NEGATIVE_P (y))
748          ? (bignum_add_unsigned (x, y, 0))
749          : (bignum_subtract_unsigned (x, y)))));
750}
751
752static bignum_type
753bignum_subtract_unsigned(bignum_type x, bignum_type y)
754{
755  int negative_p;
756  switch (bignum_compare_unsigned (x, y))
757    {
758    case bignum_comparison_equal:
759      return (BIGNUM_ZERO ());
760    case bignum_comparison_less:
761      {
762        bignum_type z = x;
763        x = y;
764        y = z;
765      }
766      negative_p = 1;
767      break;
768    case bignum_comparison_greater:
769      negative_p = 0;
770      break;
771    }
772  {
773    bignum_length_type x_length = (BIGNUM_LENGTH (x));
774    bignum_type r = (bignum_allocate (x_length, negative_p));
775    bignum_digit_type difference;
776    bignum_digit_type borrow = 0;
777    bignum_digit_type * scan_x = (BIGNUM_START_PTR (x));
778    bignum_digit_type * scan_r = (BIGNUM_START_PTR (r));
779    {
780      bignum_digit_type * scan_y = (BIGNUM_START_PTR (y));
781      bignum_digit_type * end_y = (scan_y + (BIGNUM_LENGTH (y)));
782      while (scan_y < end_y)
783        {
784          difference = (((*scan_x++) - (*scan_y++)) - borrow);
785          if (difference < 0)
786            {
787              (*scan_r++) = (difference + BIGNUM_RADIX);
788              borrow = 1;
789            }
790          else
791            {
792              (*scan_r++) = difference;
793              borrow = 0;
794            }
795        }
796    }
797    {
798      bignum_digit_type * end_x = ((BIGNUM_START_PTR (x)) + x_length);
799      if (borrow != 0)
800        while (scan_x < end_x)
801          {
802            difference = ((*scan_x++) - borrow);
803            if (difference < 0)
804              (*scan_r++) = (difference + BIGNUM_RADIX);
805            else
806              {
807                (*scan_r++) = difference;
808                borrow = 0;
809                break;
810              }
811          }
812      assert(borrow == 0);
813      while (scan_x < end_x)
814        (*scan_r++) = (*scan_x++);
815    }
816    return (bignum_trim (r));
817  }
818}
819
820
821static void
822big_times_big(C_word c, C_word self, C_word k, C_word x, C_word y)
823{
824  bignum_type bigx = big_of(x), bigy = big_of(y);
825  int neg_p = ((BIGNUM_NEGATIVE_P (bigx))
826               ? (! (BIGNUM_NEGATIVE_P (bigy)))
827               : (BIGNUM_NEGATIVE_P (bigy)));
828  /* If length 1 or 0, it should be a fixnum */
829  assert(BIGNUM_LENGTH(bigx) > 1);
830  assert(BIGNUM_LENGTH(bigy) > 1);
831  C_return_bignum(k, bignum_multiply_unsigned(bigx, bigy, neg_p));
832}
833
834static void
835fix_times_fix(C_word c, C_word self, C_word k, C_word x, C_word y)
836{
837  bignum_type bigx, bigy, result;
838  bignum_digit_type absx, absy;
839  C_word neg_p;
840
841  absx = C_unfix(x);
842  absx = absx < 0 ? -absx : absx;
843  absy = C_unfix(y);
844  absy = absy < 0 ? -absy : absy;
845  neg_p = ((x & C_INT_SIGN_BIT) ? !(y & C_INT_SIGN_BIT) : (y & C_INT_SIGN_BIT));
846
847  if (absx < BIGNUM_RADIX_ROOT) {
848     if (absx == 0 || absx == 1 || absy < BIGNUM_RADIX_ROOT) {
849       C_kontinue(k, C_fix(neg_p ? -(absx * absy) : (absx * absy)));
850     } else {
851       bigy = bignum_allocate_from_fixnum(y);
852       result = bignum_multiply_unsigned_small_factor(bigy, absx, neg_p ? 1 : 0);
853       BIGNUM_DEALLOCATE(bigy);
854       C_return_bignum(k, result);
855     }
856  } else if (absy < BIGNUM_RADIX_ROOT) {
857     if (absy == 0 || absy == 1 /*|| absx < BIGNUM_RADIX_ROOT */) {
858       C_kontinue(k, C_fix(neg_p ? -(absx * absy) : (absx * absy)));
859     } else {
860       bigx = bignum_allocate_from_fixnum(x);
861       result = bignum_multiply_unsigned_small_factor(bigx, absy, neg_p ? 1 : 0);
862       BIGNUM_DEALLOCATE(bigx);
863       C_return_bignum(k, result);
864     }
865  } else {
866    bigx = bignum_allocate_from_fixnum(x);
867    bigy = bignum_allocate_from_fixnum(y);
868    result = bignum_multiply_unsigned (bigx, bigy, neg_p ? 1 : 0);
869    BIGNUM_DEALLOCATE(bigy);
870    BIGNUM_DEALLOCATE(bigx);
871    C_return_bignum(k, result);
872  }
873}
874
875
876static void
877fix_times_big(C_word c, C_word self, C_word k, C_word x, C_word y)
878{
879  bignum_type bigx, result, bigy;
880  bignum_digit_type absx;
881  C_word neg_p;
882 
883  bigy = big_of(y);
884 
885  if (x == C_fix(0))
886    C_kontinue(k, C_fix(0));
887  else if (x == C_fix(1))
888    C_kontinue(k, y);
889  else if (x == C_fix(-1))
890    C_return_bignum(k, bignum_new_sign(bigy, !(BIGNUM_NEGATIVE_P(bigy))));
891
892  absx = C_unfix(x);
893  absx = absx < 0 ? -absx : absx; 
894  neg_p = ((x & C_INT_SIGN_BIT)
895           ? !(BIGNUM_NEGATIVE_P(bigy)) : (BIGNUM_NEGATIVE_P(bigy)));
896 
897  if (absx < BIGNUM_RADIX_ROOT) {
898     C_return_bignum(k, bignum_multiply_unsigned_small_factor(bigy, absx, neg_p));
899  } else {
900    bigx = bignum_allocate_from_fixnum(x);
901    result = bignum_multiply_unsigned (bigx, bigy, neg_p);
902    BIGNUM_DEALLOCATE(bigx);
903    C_return_bignum(k, result);
904  }
905}
906
907/* Multiplication
908   Maximum value for product_low or product_high:
909        ((R * R) + (R * (R - 2)) + (R - 1))
910   Maximum value for carry: ((R * (R - 1)) + (R - 1))
911        where R == BIGNUM_RADIX_ROOT */
912
913static bignum_type
914bignum_multiply_unsigned(bignum_type x, bignum_type y, int negative_p)
915{
916  if ((BIGNUM_LENGTH (y)) > (BIGNUM_LENGTH (x)))
917    {
918      bignum_type z = x;
919      x = y;
920      y = z;
921    }
922  {
923    bignum_digit_type carry;
924    bignum_digit_type y_digit_low;
925    bignum_digit_type y_digit_high;
926    bignum_digit_type x_digit_low;
927    bignum_digit_type x_digit_high;
928    bignum_digit_type product_low;
929    bignum_digit_type * scan_r;
930    bignum_digit_type * scan_y;
931    bignum_length_type x_length = (BIGNUM_LENGTH (x));
932    bignum_length_type y_length = (BIGNUM_LENGTH (y));
933    bignum_type r = (bignum_allocate_zeroed ((x_length + y_length), negative_p));
934    bignum_digit_type * scan_x = (BIGNUM_START_PTR (x));
935    bignum_digit_type * end_x = (scan_x + x_length);
936    bignum_digit_type * start_y = (BIGNUM_START_PTR (y));
937    bignum_digit_type * end_y = (start_y + y_length);
938    bignum_digit_type * start_r = (BIGNUM_START_PTR (r));
939#define x_digit x_digit_high
940#define y_digit y_digit_high
941#define product_high carry
942    while (scan_x < end_x)
943      {
944        x_digit = (*scan_x++);
945        x_digit_low = (HD_LOW (x_digit));
946        x_digit_high = (HD_HIGH (x_digit));
947        carry = 0;
948        scan_y = start_y;
949        scan_r = (start_r++);
950        while (scan_y < end_y)
951          {
952            y_digit = (*scan_y++);
953            y_digit_low = (HD_LOW (y_digit));
954            y_digit_high = (HD_HIGH (y_digit));
955            product_low =
956              ((*scan_r) +
957               (x_digit_low * y_digit_low) +
958               (HD_LOW (carry)));
959            product_high =
960              ((x_digit_high * y_digit_low) +
961               (x_digit_low * y_digit_high) +
962               (HD_HIGH (product_low)) +
963               (HD_HIGH (carry)));
964            (*scan_r++) =
965              (HD_CONS ((HD_LOW (product_high)), (HD_LOW (product_low))));
966            carry =
967              ((x_digit_high * y_digit_high) +
968               (HD_HIGH (product_high)));
969          }
970        (*scan_r) += carry;
971      }
972    return (bignum_trim (r));
973#undef x_digit
974#undef y_digit
975#undef product_high
976  }
977}
978
979static bignum_type
980bignum_multiply_unsigned_small_factor(bignum_type x, bignum_digit_type y,
981                                      int negative_p)
982{
983  bignum_length_type length_x = (BIGNUM_LENGTH (x));
984  bignum_type p = (bignum_allocate ((length_x + 1), negative_p));
985  bignum_destructive_copy (x, p);
986  (BIGNUM_REF (p, length_x)) = 0;
987  bignum_destructive_scale_up (p, y);
988  return (bignum_trim (p));
989}
990
991static void
992bignum_destructive_scale_up(bignum_type bignum, bignum_digit_type factor)
993{
994  bignum_digit_type carry = 0;
995  bignum_digit_type * scan = (BIGNUM_START_PTR (bignum));
996  bignum_digit_type two_digits;
997  bignum_digit_type product_low;
998#define product_high carry
999  bignum_digit_type * end = (scan + (BIGNUM_LENGTH (bignum)));
1000  assert((factor > 1) && (factor < BIGNUM_RADIX_ROOT));
1001  while (scan < end)
1002    {
1003      two_digits = (*scan);
1004      product_low = ((factor * (HD_LOW (two_digits))) + (HD_LOW (carry)));
1005      product_high =
1006        ((factor * (HD_HIGH (two_digits))) +
1007         (HD_HIGH (product_low)) +
1008         (HD_HIGH (carry)));
1009      (*scan++) = (HD_CONS ((HD_LOW (product_high)), (HD_LOW (product_low))));
1010      carry = (HD_HIGH (product_high));
1011    }
1012  /* A carry here would be an overflow, i.e. it would not fit.
1013     Hopefully the callers allocate enough space that this will
1014     never happen.
1015   */
1016  assert(carry == 0);
1017  return;
1018#undef product_high
1019}
1020
1021static void
1022bignum_destructive_add(bignum_type bignum, bignum_digit_type n)
1023{
1024  bignum_digit_type * scan = (BIGNUM_START_PTR (bignum));
1025  bignum_digit_type digit;
1026  digit = ((*scan) + n);
1027  if (digit < BIGNUM_RADIX)
1028    {
1029      (*scan) = digit;
1030      return;
1031    }
1032  (*scan++) = (digit - BIGNUM_RADIX);
1033  while (1)
1034    {
1035      digit = ((*scan) + 1);
1036      if (digit < BIGNUM_RADIX)
1037        {
1038          (*scan) = digit;
1039          return;
1040        }
1041      (*scan++) = (digit - BIGNUM_RADIX);
1042    }
1043}
1044
1045/* Division */
1046
1047/* For help understanding this algorithm, see:
1048   Knuth, Donald E., "The Art of Computer Programming",
1049   volume 2, "Seminumerical Algorithms"
1050   section 4.3.1, "Multiple-Precision Arithmetic". */
1051
1052static void
1053bignum_divide_unsigned_large_denominator(bignum_type numerator,
1054                                         bignum_type denominator,
1055                                         bignum_type * quotient,
1056                                         bignum_type * remainder,
1057                                         int q_negative_p,
1058                                         int r_negative_p)
1059{
1060  bignum_length_type length_n = ((BIGNUM_LENGTH (numerator)) + 1);
1061  bignum_length_type length_d = (BIGNUM_LENGTH (denominator));
1062  bignum_type q =
1063    ((quotient != ((bignum_type *) 0))
1064     ? (bignum_allocate ((length_n - length_d), q_negative_p))
1065     : BIGNUM_OUT_OF_BAND);
1066  bignum_type u = (bignum_allocate (length_n, r_negative_p));
1067  int shift = 0;
1068  assert(length_d > 1);
1069  {
1070    bignum_digit_type v1 = (BIGNUM_REF ((denominator), (length_d - 1)));
1071    while (v1 < (BIGNUM_RADIX / 2))
1072      {
1073        v1 <<= 1;
1074        shift += 1;
1075      }
1076  }
1077  if (shift == 0)
1078    {
1079      bignum_destructive_copy (numerator, u);
1080      (BIGNUM_REF (u, (length_n - 1))) = 0;
1081      bignum_divide_unsigned_normalized (u, denominator, q);
1082      if (remainder != ((bignum_type *) 0))
1083        (*remainder) = (bignum_trim (u));
1084      else
1085        BIGNUM_DEALLOCATE (u);
1086    }
1087  else
1088    {
1089      bignum_type v = (bignum_allocate (length_d, 0));
1090      bignum_destructive_normalization (numerator, u, shift);
1091      bignum_destructive_normalization (denominator, v, shift);
1092      bignum_divide_unsigned_normalized (u, v, q);
1093      BIGNUM_DEALLOCATE (v);
1094      if (remainder != ((bignum_type *) 0))
1095        (*remainder) = bignum_destructive_unnormalization (u, shift);
1096      else
1097        BIGNUM_DEALLOCATE(u);
1098    }
1099  if (quotient != ((bignum_type *) 0))
1100    (*quotient) = (bignum_trim (q));
1101  return;
1102}
1103
1104static void
1105bignum_divide_unsigned_normalized(bignum_type u, bignum_type v, bignum_type q)
1106{
1107  bignum_length_type u_length = (BIGNUM_LENGTH (u));
1108  bignum_length_type v_length = (BIGNUM_LENGTH (v));
1109  bignum_digit_type * u_start = (BIGNUM_START_PTR (u));
1110  bignum_digit_type * u_scan = (u_start + u_length);
1111  bignum_digit_type * u_scan_limit = (u_start + v_length);
1112  bignum_digit_type * u_scan_start = (u_scan - v_length);
1113  bignum_digit_type * v_start = (BIGNUM_START_PTR (v));
1114  bignum_digit_type * v_end = (v_start + v_length);
1115  bignum_digit_type * q_scan;
1116  bignum_digit_type v1 = (v_end[-1]);
1117  bignum_digit_type v2 = (v_end[-2]);
1118  bignum_digit_type ph; /* high half of double-digit product */
1119  bignum_digit_type pl; /* low half of double-digit product */
1120  bignum_digit_type guess;
1121  bignum_digit_type gh; /* high half-digit of guess */
1122  bignum_digit_type ch; /* high half of double-digit comparand */
1123  bignum_digit_type v2l = (HD_LOW (v2));
1124  bignum_digit_type v2h = (HD_HIGH (v2));
1125  bignum_digit_type cl; /* low half of double-digit comparand */
1126#define gl ph                   /* low half-digit of guess */
1127#define uj pl
1128#define qj ph
1129  bignum_digit_type gm;         /* memory loc for reference parameter */
1130  if (q != BIGNUM_OUT_OF_BAND)
1131    q_scan = ((BIGNUM_START_PTR (q)) + (BIGNUM_LENGTH (q)));
1132  while (u_scan_limit < u_scan)
1133    {
1134      uj = (*--u_scan);
1135      if (uj != v1)
1136        {
1137          /* comparand =
1138             (((((uj * BIGNUM_RADIX) + uj1) % v1) * BIGNUM_RADIX) + uj2);
1139             guess = (((uj * BIGNUM_RADIX) + uj1) / v1); */
1140          cl = (u_scan[-2]);
1141          ch = (bignum_digit_divide (uj, (u_scan[-1]), v1, (&gm)));
1142          guess = gm;
1143        }
1144      else
1145        {
1146          cl = (u_scan[-2]);
1147          ch = ((u_scan[-1]) + v1);
1148          guess = (BIGNUM_RADIX - 1);
1149        }
1150      while (1)
1151        {
1152          /* product = (guess * v2); */
1153          gl = (HD_LOW (guess));
1154          gh = (HD_HIGH (guess));
1155          pl = (v2l * gl);
1156          ph = ((v2l * gh) + (v2h * gl) + (HD_HIGH (pl)));
1157          pl = (HD_CONS ((HD_LOW (ph)), (HD_LOW (pl))));
1158          ph = ((v2h * gh) + (HD_HIGH (ph)));
1159          /* if (comparand >= product) */
1160          if ((ch > ph) || ((ch == ph) && (cl >= pl)))
1161            break;
1162          guess -= 1;
1163          /* comparand += (v1 << BIGNUM_DIGIT_LENGTH) */
1164          ch += v1;
1165          /* if (comparand >= (BIGNUM_RADIX * BIGNUM_RADIX)) */
1166          if (ch >= BIGNUM_RADIX)
1167            break;
1168        }
1169      qj = (bignum_divide_subtract (v_start, v_end, guess, (--u_scan_start)));
1170      if (q != BIGNUM_OUT_OF_BAND)
1171        (*--q_scan) = qj;
1172    }
1173  return;
1174#undef gl
1175#undef uj
1176#undef qj
1177}
1178
1179static bignum_digit_type
1180bignum_divide_subtract(bignum_digit_type * v_start,
1181                       bignum_digit_type * v_end,
1182                       bignum_digit_type guess,
1183                       bignum_digit_type * u_start)
1184{
1185  bignum_digit_type * v_scan = v_start;
1186  bignum_digit_type * u_scan = u_start;
1187  bignum_digit_type carry = 0;
1188  if (guess == 0) return (0);
1189  {
1190    bignum_digit_type gl = (HD_LOW (guess));
1191    bignum_digit_type gh = (HD_HIGH (guess));
1192    bignum_digit_type v;
1193    bignum_digit_type pl;
1194    bignum_digit_type vl;
1195#define vh v
1196#define ph carry
1197#define diff pl
1198    while (v_scan < v_end)
1199      {
1200        v = (*v_scan++);
1201        vl = (HD_LOW (v));
1202        vh = (HD_HIGH (v));
1203        pl = ((vl * gl) + (HD_LOW (carry)));
1204        ph = ((vl * gh) + (vh * gl) + (HD_HIGH (pl)) + (HD_HIGH (carry)));
1205        diff = ((*u_scan) - (HD_CONS ((HD_LOW (ph)), (HD_LOW (pl)))));
1206        if (diff < 0)
1207          {
1208            (*u_scan++) = (diff + BIGNUM_RADIX);
1209            carry = ((vh * gh) + (HD_HIGH (ph)) + 1);
1210          }
1211        else
1212          {
1213            (*u_scan++) = diff;
1214            carry = ((vh * gh) + (HD_HIGH (ph)));
1215          }
1216      }
1217    if (carry == 0)
1218      return (guess);
1219    diff = ((*u_scan) - carry);
1220    if (diff < 0)
1221      (*u_scan) = (diff + BIGNUM_RADIX);
1222    else
1223      {
1224        (*u_scan) = diff;
1225        return (guess);
1226      }
1227#undef vh
1228#undef ph
1229#undef diff
1230  }
1231  /* Subtraction generated carry, implying guess is one too large.
1232     Add v back in to bring it back down. */
1233  v_scan = v_start;
1234  u_scan = u_start;
1235  carry = 0;
1236  while (v_scan < v_end)
1237    {
1238      bignum_digit_type sum = ((*v_scan++) + (*u_scan) + carry);
1239      if (sum < BIGNUM_RADIX)
1240        {
1241          (*u_scan++) = sum;
1242          carry = 0;
1243        }
1244      else
1245        {
1246          (*u_scan++) = (sum - BIGNUM_RADIX);
1247          carry = 1;
1248        }
1249    }
1250  if (carry == 1)
1251    {
1252      bignum_digit_type sum = ((*u_scan) + carry);
1253      (*u_scan) = ((sum < BIGNUM_RADIX) ? sum : (sum - BIGNUM_RADIX));
1254    }
1255  return (guess - 1);
1256}
1257
1258static void
1259bignum_divide_unsigned_medium_denominator(bignum_type numerator,
1260                                          bignum_digit_type denominator,
1261                                          bignum_type * quotient,
1262                                          bignum_type * remainder,
1263                                          int q_negative_p,
1264                                          int r_negative_p)
1265{
1266  bignum_length_type length_n = (BIGNUM_LENGTH (numerator));
1267  bignum_length_type length_q;
1268  bignum_type q;
1269  int shift = 0;
1270  /* Because `bignum_digit_divide' requires a normalized denominator. */
1271  while (denominator < (BIGNUM_RADIX / 2))
1272    {
1273      denominator <<= 1;
1274      shift += 1;
1275    }
1276  if (shift == 0)
1277    {
1278      length_q = length_n;
1279      q = (bignum_allocate (length_q, q_negative_p));
1280      bignum_destructive_copy (numerator, q);
1281    }
1282  else
1283    {
1284      length_q = (length_n + 1);
1285      q = (bignum_allocate (length_q, q_negative_p));
1286      bignum_destructive_normalization (numerator, q, shift);
1287    }
1288  {
1289    bignum_digit_type r = 0;
1290    bignum_digit_type * start = (BIGNUM_START_PTR (q));
1291    bignum_digit_type * scan = (start + length_q);
1292    bignum_digit_type qj;
1293    if (quotient != ((bignum_type *) 0))
1294      {
1295        while (start < scan)
1296          {
1297            r = (bignum_digit_divide (r, (*--scan), denominator, (&qj)));
1298            (*scan) = qj;
1299          }
1300        (*quotient) = (bignum_trim (q));
1301      }
1302    else
1303      {
1304        while (start < scan)
1305          r = (bignum_digit_divide (r, (*--scan), denominator, (&qj)));
1306        BIGNUM_DEALLOCATE (q);
1307      }
1308    if (remainder != ((bignum_type *) 0))
1309      {
1310        if (shift != 0)
1311          r >>= shift;
1312        (*remainder) = (bignum_digit_to_bignum (r, r_negative_p));
1313      }
1314  }
1315  return;
1316}
1317
1318static void
1319bignum_destructive_normalization(bignum_type source, bignum_type target,
1320                                 int shift_left)
1321{
1322  bignum_digit_type digit;
1323  bignum_digit_type * scan_source = (BIGNUM_START_PTR (source));
1324  bignum_digit_type carry = 0;
1325  bignum_digit_type * scan_target = (BIGNUM_START_PTR (target));
1326  bignum_digit_type * end_source = (scan_source + (BIGNUM_LENGTH (source)));
1327  bignum_digit_type * end_target = (scan_target + (BIGNUM_LENGTH (target)));
1328  int shift_right = (BIGNUM_DIGIT_LENGTH - shift_left);
1329  bignum_digit_type mask = ((1L << shift_right) - 1);
1330  while (scan_source < end_source)
1331    {
1332      digit = (*scan_source++);
1333      (*scan_target++) = (((digit & mask) << shift_left) | carry);
1334      carry = (digit >> shift_right);
1335    }
1336  if (scan_target < end_target)
1337    (*scan_target) = carry;
1338  else
1339    assert(carry == 0);
1340  return;
1341}
1342
1343/* This will also trim the number if necessary */
1344static bignum_type
1345bignum_destructive_unnormalization(bignum_type bignum, int shift_right)
1346{
1347  bignum_length_type length = BIGNUM_LENGTH(bignum);
1348  bignum_digit_type * start = (BIGNUM_START_PTR (bignum));
1349  bignum_digit_type * scan = (start + length);
1350  bignum_digit_type digit;
1351  bignum_digit_type carry = 0;
1352  int shift_left = (BIGNUM_DIGIT_LENGTH - shift_right);
1353  bignum_digit_type mask = ((1L << shift_right) - 1);
1354 
1355  while (!(*--scan)) {
1356    if (start == scan) { /* Don't bother. */
1357      BIGNUM_SET_HEADER (bignum, 0, 0);
1358      BIGNUM_REDUCE_LENGTH (bignum, bignum, 0);
1359      return bignum;
1360    }
1361    --length;
1362  }
1363
1364  digit = (*scan);
1365  (*scan) = (digit >> shift_right);
1366  length -= (*scan == 0); /* Add 1 or 0 */
1367  carry = ((digit & mask) << shift_left);
1368 
1369  while (start < scan)
1370    {
1371      digit = (*--scan);
1372      (*scan) = ((digit >> shift_right) | carry);
1373      carry = ((digit & mask) << shift_left);
1374    }
1375  assert(carry == 0);
1376  assert(BIGNUM_LENGTH(bignum) != length);
1377  assert(length != 1 || BIGNUM_REF(bignum, 0) != 0);
1378  BIGNUM_SET_HEADER
1379    (bignum, length, (BIGNUM_NEGATIVE_P (bignum)));
1380  BIGNUM_REDUCE_LENGTH (bignum, bignum, length);
1381
1382  return bignum;
1383}
1384
1385/* This is a reduced version of the division algorithm, applied to the
1386   case of dividing two bignum digits by one bignum digit.  It is
1387   assumed that the numerator, denominator are normalized. */
1388
1389#define BDD_STEP(qn, j)                                                 \
1390{                                                                       \
1391  uj = (u[j]);                                                          \
1392  if (uj != v1)                                                         \
1393    {                                                                   \
1394      uj_uj1 = (HD_CONS (uj, (u[j + 1])));                              \
1395      guess = (uj_uj1 / v1);                                            \
1396      comparand = (HD_CONS ((uj_uj1 % v1), (u[j + 2])));                \
1397    }                                                                   \
1398  else                                                                  \
1399    {                                                                   \
1400      guess = (BIGNUM_RADIX_ROOT - 1);                                  \
1401      comparand = (HD_CONS (((u[j + 1]) + v1), (u[j + 2])));            \
1402    }                                                                   \
1403  while ((guess * v2) > comparand)                                      \
1404    {                                                                   \
1405      guess -= 1;                                                       \
1406      comparand += (v1 << BIGNUM_HALF_DIGIT_LENGTH);                    \
1407      if (comparand >= BIGNUM_RADIX)                                    \
1408        break;                                                          \
1409    }                                                                   \
1410  qn = (bignum_digit_divide_subtract (v1, v2, guess, (&u[j])));         \
1411}
1412
1413static bignum_digit_type
1414bignum_digit_divide(bignum_digit_type uh, bignum_digit_type ul,
1415                    bignum_digit_type v,
1416                    bignum_digit_type * q) /* return value */
1417{
1418  bignum_digit_type guess;
1419  bignum_digit_type comparand;
1420  bignum_digit_type v1;
1421  bignum_digit_type v2;
1422  bignum_digit_type uj;
1423  bignum_digit_type uj_uj1;
1424  bignum_digit_type q1;
1425  bignum_digit_type q2;
1426  bignum_digit_type u [4];
1427  if (uh == 0)
1428    {
1429      if (ul < v)
1430        {
1431          (*q) = 0;
1432          return (ul);
1433        }
1434      else if (ul == v)
1435        {
1436          (*q) = 1;
1437          return (0);
1438        }
1439    }
1440  (u[0]) = (HD_HIGH (uh));
1441  (u[1]) = (HD_LOW (uh));
1442  (u[2]) = (HD_HIGH (ul));
1443  (u[3]) = (HD_LOW (ul));
1444  v1 = (HD_HIGH (v));
1445  v2 = (HD_LOW (v));
1446  BDD_STEP (q1, 0);
1447  BDD_STEP (q2, 1);
1448  (*q) = (HD_CONS (q1, q2));
1449  return (HD_CONS ((u[2]), (u[3])));
1450}
1451
1452#undef BDD_STEP
1453
1454#define BDDS_MULSUB(vn, un, carry_in)                                   \
1455{                                                                       \
1456  product = ((vn * guess) + carry_in);                                  \
1457  diff = (un - (HD_LOW (product)));                                     \
1458  if (diff < 0)                                                         \
1459    {                                                                   \
1460      un = (diff + BIGNUM_RADIX_ROOT);                                  \
1461      carry = ((HD_HIGH (product)) + 1);                                \
1462    }                                                                   \
1463  else                                                                  \
1464    {                                                                   \
1465      un = diff;                                                        \
1466      carry = (HD_HIGH (product));                                      \
1467    }                                                                   \
1468}
1469
1470#define BDDS_ADD(vn, un, carry_in)                                      \
1471{                                                                       \
1472  sum = (vn + un + carry_in);                                           \
1473  if (sum < BIGNUM_RADIX_ROOT)                                          \
1474    {                                                                   \
1475      un = sum;                                                         \
1476      carry = 0;                                                        \
1477    }                                                                   \
1478  else                                                                  \
1479    {                                                                   \
1480      un = (sum - BIGNUM_RADIX_ROOT);                                   \
1481      carry = 1;                                                        \
1482    }                                                                   \
1483}
1484
1485static bignum_digit_type
1486bignum_digit_divide_subtract(bignum_digit_type v1, bignum_digit_type v2,
1487                             bignum_digit_type guess, bignum_digit_type * u)
1488{
1489  {
1490    bignum_digit_type product;
1491    bignum_digit_type diff;
1492    bignum_digit_type carry;
1493    BDDS_MULSUB (v2, (u[2]), 0);
1494    BDDS_MULSUB (v1, (u[1]), carry);
1495    if (carry == 0)
1496      return (guess);
1497    diff = ((u[0]) - carry);
1498    if (diff < 0)
1499      (u[0]) = (diff + BIGNUM_RADIX);
1500    else
1501      {
1502        (u[0]) = diff;
1503        return (guess);
1504      }
1505  }
1506  {
1507    bignum_digit_type sum;
1508    bignum_digit_type carry;
1509    BDDS_ADD(v2, (u[2]), 0);
1510    BDDS_ADD(v1, (u[1]), carry);
1511    if (carry == 1)
1512      (u[0]) += 1;
1513  }
1514  return (guess - 1);
1515}
1516
1517#undef BDDS_MULSUB
1518#undef BDDS_ADD
1519
1520static void
1521bignum_divide_unsigned_small_denominator(bignum_type numerator,
1522                                         bignum_digit_type denominator,
1523                                         bignum_type * quotient,
1524                                         bignum_type * remainder,
1525                                         int q_negative_p,  int r_negative_p)
1526{
1527  bignum_type q = (bignum_new_sign (numerator, q_negative_p));
1528  bignum_digit_type r = (bignum_destructive_scale_down (q, denominator));
1529  (*quotient) = (bignum_trim (q));
1530  if (remainder != ((bignum_type *) 0))
1531    (*remainder) = (bignum_digit_to_bignum (r, r_negative_p));
1532  return;
1533}
1534
1535/* Given (denominator > 1), it is fairly easy to show that
1536   (quotient_high < BIGNUM_RADIX_ROOT), after which it is easy to see
1537   that all digits are < BIGNUM_RADIX. */
1538
1539static bignum_digit_type
1540bignum_destructive_scale_down(bignum_type bignum, bignum_digit_type denominator)
1541{
1542  bignum_digit_type numerator;
1543  bignum_digit_type remainder = 0;
1544  bignum_digit_type two_digits;
1545#define quotient_high remainder
1546  bignum_digit_type * start = (BIGNUM_START_PTR (bignum));
1547  bignum_digit_type * scan = (start + (BIGNUM_LENGTH (bignum)));
1548  assert((denominator > 1) && (denominator < BIGNUM_RADIX_ROOT));
1549  while (start < scan)
1550    {
1551      two_digits = (*--scan);
1552      numerator = (HD_CONS (remainder, (HD_HIGH (two_digits))));
1553      quotient_high = (numerator / denominator);
1554      numerator = (HD_CONS ((numerator % denominator), (HD_LOW (two_digits))));
1555      (*scan) = (HD_CONS (quotient_high, (numerator / denominator)));
1556      remainder = (numerator % denominator);
1557    }
1558  return (remainder);
1559#undef quotient_high
1560}
1561
1562static bignum_type
1563bignum_remainder_unsigned_small_denominator(bignum_type n, bignum_digit_type d,
1564                                            int negative_p)
1565{
1566  bignum_digit_type two_digits;
1567  bignum_digit_type * start = (BIGNUM_START_PTR (n));
1568  bignum_digit_type * scan = (start + (BIGNUM_LENGTH (n)));
1569  bignum_digit_type r = 0;
1570  assert((d > 1) && (d < BIGNUM_RADIX_ROOT));
1571  while (start < scan)
1572    {
1573      two_digits = (*--scan);
1574      r =
1575        ((HD_CONS (((HD_CONS (r, (HD_HIGH (two_digits)))) % d),
1576                   (HD_LOW (two_digits))))
1577         % d);
1578    }
1579  return (bignum_digit_to_bignum (r, negative_p));
1580}
1581
1582static C_word
1583big_comp_big(C_word x, C_word y)
1584{
1585  bignum_type bigx = big_of(x), bigy = big_of(y);
1586  return C_fix((BIGNUM_NEGATIVE_P (bigx))
1587               ? ((BIGNUM_NEGATIVE_P (bigy))
1588                  ? (bignum_compare_unsigned (bigy, bigx))
1589                  : (bignum_comparison_less))
1590               : ((BIGNUM_NEGATIVE_P (bigy))
1591                  ? (bignum_comparison_greater)
1592                  : (bignum_compare_unsigned (bigx, bigy))));
1593}
1594
1595static enum bignum_comparison
1596bignum_compare_unsigned(bignum_type x, bignum_type y)
1597{
1598  bignum_length_type x_length;
1599  bignum_length_type y_length;
1600  if (x == y) /* Objects are the same? */
1601    return (bignum_comparison_equal);
1602
1603  x_length = (BIGNUM_LENGTH (x));
1604  y_length = (BIGNUM_LENGTH (y));
1605  if (x_length < y_length)
1606    return (bignum_comparison_less);
1607  if (x_length > y_length)
1608    return (bignum_comparison_greater);
1609  {
1610    bignum_digit_type * start_x = (BIGNUM_START_PTR (x));
1611    bignum_digit_type * scan_x = (start_x + x_length);
1612    bignum_digit_type * scan_y = ((BIGNUM_START_PTR (y)) + y_length);
1613    while (start_x < scan_x)
1614      {
1615        bignum_digit_type digit_x = (*--scan_x);
1616        bignum_digit_type digit_y = (*--scan_y);
1617        if (digit_x < digit_y)
1618          return (bignum_comparison_less);
1619        if (digit_x > digit_y)
1620          return (bignum_comparison_greater);
1621      }
1622  }
1623  return (bignum_comparison_equal);
1624}
1625
1626static void
1627big_abs(C_word c, C_word self, C_word k, C_word big)
1628{
1629  if (!BIGNUM_NEGATIVE_P(big_of(big)))
1630    C_kontinue(k, big);
1631  else
1632    C_return_bignum(k, bignum_new_sign(big_of(big), 0));
1633}
1634
1635static void
1636big_quotient_fix(C_word c, C_word self, C_word k, C_word x, C_word y)
1637{
1638  bignum_type bigx = big_of(x);
1639  y = C_unfix(y);
1640
1641  if (y == 1)
1642    C_kontinue(k, x);
1643  else if (y == -1)
1644    C_return_bignum(k, bignum_new_sign(bigx, !(BIGNUM_NEGATIVE_P(bigx))));
1645
1646  /* Too bad, we really need to do some work... */
1647  {
1648    int neg_p = (y < 0) ? !(BIGNUM_NEGATIVE_P(bigx)) : BIGNUM_NEGATIVE_P(bigx);
1649    bignum_digit_type abs_y = (y < 0) ? -y : y;
1650    bignum_type quotient;
1651   
1652    if (y == C_MOST_NEGATIVE_FIXNUM) {
1653      if (!BIGNUM_NEGATIVE_P(bigx) && BIGNUM_LENGTH(bigx) == 1
1654          && BIGNUM_REF(bigx, 1) == 1 && BIGNUM_REF(bigx, 0) == 0) {
1655        /*
1656         * Very very special case:
1657         * quotient(MOST_NEGATIVE_FIXNUM, -(MOST_NEGATIVE_FIXNUM)) => -1
1658         */
1659        C_kontinue(k, C_fix(-1));
1660      } else {
1661        /* This is the only case we need to go allocate a bignum for */
1662        bignum_type bigy =
1663          bignum_allocate_from_fixnum(C_fix(C_MOST_NEGATIVE_FIXNUM));
1664
1665        bignum_divide_unsigned_large_denominator
1666          (bigx, bigy, (&quotient), ((bignum_type *) 0), neg_p, 0);
1667        BIGNUM_DEALLOCATE(bigy);
1668        C_return_bignum(k, quotient);
1669      }
1670    } else if (abs_y < BIGNUM_RADIX_ROOT) {
1671      bignum_divide_unsigned_small_denominator
1672        (bigx, abs_y, (&quotient), ((bignum_type *) 0), neg_p, 0);
1673      C_return_bignum(k, quotient);
1674    } else {
1675      bignum_divide_unsigned_medium_denominator
1676        (bigx, abs_y, (&quotient), ((bignum_type *) 0), neg_p, 0);
1677      C_return_bignum(k, quotient);
1678    }
1679  }
1680}
1681
1682static void
1683big_quotient_big(C_word c, C_word self, C_word k, C_word x, C_word y)
1684{
1685  bignum_type numerator = big_of(x);
1686  bignum_type denominator = big_of(y);
1687  int neg_p = ((BIGNUM_NEGATIVE_P (denominator))
1688               ? (! (BIGNUM_NEGATIVE_P (numerator)))
1689               : (BIGNUM_NEGATIVE_P (numerator)));
1690
1691  switch (bignum_compare_unsigned (numerator, denominator))
1692    {
1693    case bignum_comparison_equal:
1694      C_kontinue(k, neg_p ? C_fix(-1) : C_fix(1));
1695    case bignum_comparison_less:
1696      C_kontinue(k, C_fix(0));
1697    case bignum_comparison_greater:
1698    default:                                  /* to appease gcc -Wall */
1699      {
1700        bignum_type quotient;
1701        assert(BIGNUM_LENGTH(denominator) > 1);
1702        bignum_divide_unsigned_large_denominator
1703          (numerator, denominator, (&quotient), ((bignum_type *) 0), neg_p, 0);
1704        C_return_bignum(k, quotient);
1705      }
1706    }
1707}
1708
1709static void
1710big_remainder_fix(C_word c, C_word self, C_word k, C_word x, C_word y)
1711{
1712   bignum_type remainder, bigy, bigx;
1713   int x_neg, y_neg;
1714   C_word abs_y;
1715
1716   switch (y) {
1717   /* case 0:  SHOULD NOT HAPPEN (checked in Scheme)
1718     C_kontinue(k, C_SCHEME_FALSE); */
1719   case C_fix(1):
1720   case C_fix(-1):
1721     C_kontinue(k, C_fix(0));
1722   case C_fix(C_MOST_NEGATIVE_FIXNUM):
1723     bigx = big_of(x);
1724     if (BIGNUM_LENGTH(bigx) == 2 && BIGNUM_REF(bigx, 0) == 0 && BIGNUM_REF(bigx, 0))
1725       C_kontinue(k, C_fix(0));
1726     /* Don't handle 0 <= length(bigx) <= 1 since then it should be a fixnum */
1727     assert(BIGNUM_LENGTH(bigx) >= 2);
1728     
1729     bigy = bignum_allocate_from_fixnum(y);
1730     bignum_divide_unsigned_large_denominator
1731      (bigx, bigy, (bignum_type *)0, &remainder, 0, BIGNUM_NEGATIVE_P(bigx));
1732     BIGNUM_DEALLOCATE(bigy);
1733     
1734     C_return_bignum(k, remainder);
1735   default:
1736     bigx = big_of(x);
1737     y = C_unfix(y);
1738     y_neg = (y < 0);
1739     abs_y = y_neg ? -y : y;
1740     x_neg = BIGNUM_NEGATIVE_P(bigx);
1741     
1742     if (abs_y < BIGNUM_RADIX_ROOT)
1743       remainder =
1744         bignum_remainder_unsigned_small_denominator(bigx, abs_y, x_neg);
1745     else
1746       bignum_divide_unsigned_medium_denominator(bigx, abs_y,
1747                                                 (bignum_type *) 0, &remainder,
1748                                                 x_neg, x_neg);
1749     C_return_bignum(k, remainder);
1750   }
1751}
1752
1753static void
1754big_remainder_big(C_word c, C_word self, C_word k, C_word x, C_word y)
1755{
1756  bignum_type numerator = big_of(x), denominator = big_of(y);
1757 
1758  switch (bignum_compare_unsigned (numerator, denominator))
1759    {
1760      case bignum_comparison_equal:
1761        C_kontinue(k, C_fix(0));
1762      case bignum_comparison_less:
1763        C_kontinue(k, x);
1764      case bignum_comparison_greater:
1765      default:                                  /* to appease gcc -Wall */
1766        {
1767          bignum_type remainder;
1768          bignum_divide_unsigned_large_denominator
1769            (numerator, denominator,
1770             ((bignum_type *) 0), (&remainder),
1771             0, BIGNUM_NEGATIVE_P(numerator));
1772          C_return_bignum(k, remainder);
1773        }
1774    }
1775}
1776
1777static void
1778flo_to_int(C_word c, C_word self, C_word k, C_word x)
1779{
1780  C_return_bignum(k, double_to_bignum(C_flonum_magnitude(x)));
1781}
1782
1783#define DTB_WRITE_DIGIT(factor)                                         \
1784{                                                                       \
1785  significand *= (factor);                                              \
1786  digit = ((bignum_digit_type) significand);                            \
1787  (*--scan) = digit;                                                    \
1788  significand -= ((double) digit);                                      \
1789}
1790
1791static bignum_type
1792double_to_bignum(double x)
1793{
1794  int exponent;
1795  double significand = (frexp (x, (&exponent)));
1796  if (exponent <= 0) return (BIGNUM_ZERO ());
1797  if (exponent == 1) return (BIGNUM_ONE (x < 0));
1798  if (significand < 0) significand = (-significand);
1799  {
1800    bignum_length_type length = (BIGNUM_BITS_TO_DIGITS (exponent));
1801    bignum_type result = (bignum_allocate (length, (x < 0)));
1802    bignum_digit_type * start = (BIGNUM_START_PTR (result));
1803    bignum_digit_type * scan = (start + length);
1804    bignum_digit_type digit;
1805    int odd_bits = (exponent % BIGNUM_DIGIT_LENGTH);
1806    if (odd_bits > 0)
1807      DTB_WRITE_DIGIT (1L << odd_bits);
1808    while (start < scan)
1809      {
1810        if (significand == 0)
1811          {
1812            while (start < scan)
1813              (*--scan) = 0;
1814            break;
1815          }
1816        DTB_WRITE_DIGIT (BIGNUM_RADIX);
1817      }
1818    return (result);
1819  }
1820}
1821
1822#undef DTB_WRITE_DIGIT
1823
1824static void 
1825int_bitwise_int(C_word c, C_word self, C_word k, C_word op, C_word x, C_word y)
1826{
1827  bignum_type bigx, bigy, result;
1828 
1829  if((x & C_FIXNUM_BIT) != 0)
1830    bigx = bignum_allocate_from_fixnum(x);
1831  else
1832    bigx = big_of(x);
1833
1834  if((y & C_FIXNUM_BIT) != 0)
1835    bigy = bignum_allocate_from_fixnum(y);
1836  else
1837    bigy = big_of(y);
1838
1839  result = ((BIGNUM_NEGATIVE_P (bigx))
1840             ? (BIGNUM_NEGATIVE_P (bigy))
1841               ? bignum_negneg_bitwise_op(C_unfix(op), bigx, bigy)
1842               : bignum_posneg_bitwise_op(C_unfix(op), bigy, bigx)
1843             : (BIGNUM_NEGATIVE_P (bigy))
1844               ? bignum_posneg_bitwise_op(C_unfix(op), bigx, bigy)
1845               : bignum_pospos_bitwise_op(C_unfix(op), bigx, bigy));
1846       
1847  if((x & C_FIXNUM_BIT) != 0)
1848    BIGNUM_DEALLOCATE(bigx);
1849  if((y & C_FIXNUM_BIT) != 0)
1850    BIGNUM_DEALLOCATE(bigy);
1851
1852  C_return_bignum(k, result);
1853}
1854
1855static void 
1856int_not(C_word c, C_word self, C_word k, C_word x)
1857{
1858  bignum_type bigx, result;
1859 
1860  if((x & C_FIXNUM_BIT) != 0)
1861    bigx = bignum_allocate_from_fixnum(x);
1862  else
1863    bigx = big_of(x);
1864
1865  result = bignum_bitwise_not(bigx);
1866
1867  if((x & C_FIXNUM_BIT) != 0)
1868    BIGNUM_DEALLOCATE(bigx);
1869
1870  C_return_bignum(k, result);
1871}
1872
1873bignum_type
1874bignum_bitwise_not(bignum_type x)
1875{
1876  static C_word invbits[] = { BIGNUM_RADIX | 1, 1 }; /* bignum representing -1 */
1877  return bignum_subtract((bignum_type)invbits, x);
1878}
1879
1880/*
1881 * From Hacker's Delight by Frank Warren
1882 * based on a modified nlz() from section 5-3 (fig. 5-7)
1883 */
1884static C_word
1885ilen(C_uword x) {
1886  C_uword y;
1887  C_word n = 0;
1888
1889#ifdef C_SIXTY_FOUR
1890  y = x >> 32; if (y != 0) { n += 32; x = y; }
1891#endif
1892  y = x >> 16; if (y != 0) { n += 16; x = y; }
1893  y = x >>  8; if (y != 0) { n +=  8; x = y; }
1894  y = x >>  4; if (y != 0) { n +=  4; x = y; }
1895  y = x >>  2; if (y != 0) { n +=  2; x = y; }
1896  y = x >>  1; if (y != 0) return n + 2;
1897  return n + x;
1898}
1899
1900static void
1901int_length(C_word c, C_word self, C_word k, C_word x)
1902{
1903  if (x & C_FIXNUM_BIT) {
1904    x = C_unfix(x);
1905    C_kontinue(k, C_fix(ilen((x < 0) ? ~x : x)));
1906  } else { /* bignum */
1907    bignum_type bigx;
1908    bignum_digit_type bd;
1909    int len_1;
1910
1911    bigx = big_of(x);
1912    if (BIGNUM_NEGATIVE_P(bigx)) {
1913      /* TODO: This can probably be done faster, by juggling with borrow */
1914      bignum_type notx = bignum_bitwise_not(bigx);
1915      len_1 = BIGNUM_LENGTH(notx) - 1;
1916      bd = *(BIGNUM_START_PTR(notx) + len_1); /* Most significant digit */
1917      BIGNUM_DEALLOCATE(notx);
1918    } else {
1919      len_1 = BIGNUM_LENGTH(bigx) - 1;
1920      bd = *(BIGNUM_START_PTR(bigx) + len_1); /* Most significant digit */
1921    }
1922    C_kontinue(k, C_fix(ilen(bd) + len_1 * BIGNUM_DIGIT_LENGTH));
1923  }
1924}
1925
1926static void 
1927int_shift_fix(C_word c, C_word self, C_word k, C_word x, C_word y)
1928{
1929  bignum_type bigx, result;
1930
1931  if (y == C_fix(0)) C_kontinue(k, x); /* Done too (no shift) */
1932 
1933  /* Ensure x is a bignum */
1934  if((x & C_FIXNUM_BIT) != 0) {
1935    if (x == C_fix(0)) /* Skip everything else */
1936      C_kontinue(k, x);
1937   
1938    bigx = bignum_allocate_from_fixnum(x);
1939  } else {
1940    bigx = big_of(x);
1941  }
1942   
1943  result = bignum_arithmetic_shift(bigx, C_unfix(y));
1944  if ((x & C_FIXNUM_BIT) != 0)
1945    BIGNUM_DEALLOCATE(bigx);
1946  C_return_bignum(k, result);
1947}
1948
1949static bignum_type
1950bignum_arithmetic_shift(bignum_type arg1, C_word n)
1951{
1952  bignum_type tmp1, tmp2, result;
1953  if (BIGNUM_NEGATIVE_P(arg1) && n < 0) {
1954    tmp1 = bignum_bitwise_not(arg1);
1955    tmp2 = bignum_magnitude_ash(tmp1, n);
1956    BIGNUM_DEALLOCATE(tmp1);
1957    result = bignum_bitwise_not(tmp2);
1958    BIGNUM_DEALLOCATE(tmp2);
1959    return result;
1960  } else {
1961    return bignum_magnitude_ash(arg1, n);
1962  }
1963}
1964
1965static bignum_type
1966bignum_magnitude_ash(bignum_type arg1, C_word n)
1967{
1968  bignum_type result;
1969  bignum_digit_type *scan1;
1970  bignum_digit_type *scanr;
1971  bignum_digit_type *end;
1972
1973  C_word digit_offset,bit_offset;
1974  assert(n != 0); /* int_shift_fix checks this */
1975
1976  if (n > 0) {
1977    digit_offset = n / BIGNUM_DIGIT_LENGTH;
1978    bit_offset =   n % BIGNUM_DIGIT_LENGTH;
1979   
1980    result = bignum_allocate_zeroed (BIGNUM_LENGTH (arg1) + digit_offset + 1,
1981                                     BIGNUM_NEGATIVE_P(arg1));
1982
1983    scanr = BIGNUM_START_PTR (result) + digit_offset;
1984    scan1 = BIGNUM_START_PTR (arg1);
1985    end = scan1 + BIGNUM_LENGTH (arg1);
1986   
1987    while (scan1 < end) {
1988      *scanr = *scanr | (*scan1 & BIGNUM_DIGIT_MASK) << bit_offset;
1989      *scanr = *scanr & BIGNUM_DIGIT_MASK;
1990      scanr++;
1991      *scanr = *scan1++ >> (BIGNUM_DIGIT_LENGTH - bit_offset);
1992      *scanr = *scanr & BIGNUM_DIGIT_MASK;
1993    }
1994  }
1995  else if (n < 0
1996           && (-n >= (BIGNUM_LENGTH (arg1) * (bignum_length_type) BIGNUM_DIGIT_LENGTH)))
1997    result = BIGNUM_ZERO ();
1998
1999  else /* if (n < 0) */ {
2000    digit_offset = -n / BIGNUM_DIGIT_LENGTH;
2001    bit_offset =   -n % BIGNUM_DIGIT_LENGTH;
2002   
2003    result = bignum_allocate_zeroed (BIGNUM_LENGTH (arg1) - digit_offset,
2004                                     BIGNUM_NEGATIVE_P(arg1));
2005   
2006    scanr = BIGNUM_START_PTR (result);
2007    scan1 = BIGNUM_START_PTR (arg1) + digit_offset;
2008    end = scanr + BIGNUM_LENGTH (result) - 1;
2009   
2010    while (scanr < end) {
2011      *scanr =  (*scan1++ & BIGNUM_DIGIT_MASK) >> bit_offset ;
2012      *scanr = (*scanr | 
2013        *scan1 << (BIGNUM_DIGIT_LENGTH - bit_offset)) & BIGNUM_DIGIT_MASK;
2014      scanr++;
2015    }
2016    *scanr =  (*scan1++ & BIGNUM_DIGIT_MASK) >> bit_offset ;
2017  }
2018 
2019  return (bignum_trim (result));
2020}
2021
2022static bignum_type
2023bignum_pospos_bitwise_op(int op, bignum_type arg1, bignum_type arg2)
2024{
2025  bignum_type result;
2026  bignum_length_type max_length;
2027
2028  bignum_digit_type *scan1, *end1, digit1;
2029  bignum_digit_type *scan2, *end2, digit2;
2030  bignum_digit_type *scanr, *endr;
2031
2032  max_length =  (BIGNUM_LENGTH(arg1) > BIGNUM_LENGTH(arg2))
2033               ? BIGNUM_LENGTH(arg1) : BIGNUM_LENGTH(arg2);
2034
2035  result = bignum_allocate(max_length, 0);
2036
2037  scanr = BIGNUM_START_PTR(result);
2038  scan1 = BIGNUM_START_PTR(arg1);
2039  scan2 = BIGNUM_START_PTR(arg2);
2040  endr = scanr + max_length;
2041  end1 = scan1 + BIGNUM_LENGTH(arg1);
2042  end2 = scan2 + BIGNUM_LENGTH(arg2);
2043
2044  while (scanr < endr) {
2045    digit1 = (scan1 < end1) ? *scan1++ : 0;
2046    digit2 = (scan2 < end2) ? *scan2++ : 0;
2047    /*
2048    fprintf(stderr, "[pospos op = %d, i = %ld, d1 = %lx, d2 = %lx]\n",
2049            op, endr - scanr, digit1, digit2);
2050            */
2051    *scanr++ = (op == bignum_and_op) ? digit1 & digit2 :
2052               (op == bignum_ior_op) ? digit1 | digit2 :
2053                                       digit1 ^ digit2;
2054  }
2055  return bignum_trim(result);
2056}
2057
2058static bignum_type
2059bignum_posneg_bitwise_op(int op, bignum_type arg1, bignum_type arg2)
2060{
2061  bignum_type result;
2062  bignum_length_type max_length;
2063
2064  bignum_digit_type *scan1, *end1, digit1;
2065  bignum_digit_type *scan2, *end2, digit2, carry2;
2066  bignum_digit_type *scanr, *endr;
2067
2068  char neg_p = op == bignum_ior_op || op == bignum_xor_op;
2069
2070  max_length =  (BIGNUM_LENGTH(arg1) > BIGNUM_LENGTH(arg2) + 1)
2071               ? BIGNUM_LENGTH(arg1) : BIGNUM_LENGTH(arg2) + 1;
2072
2073  result = bignum_allocate(max_length, neg_p);
2074
2075  scanr = BIGNUM_START_PTR(result);
2076  scan1 = BIGNUM_START_PTR(arg1);
2077  scan2 = BIGNUM_START_PTR(arg2);
2078  endr = scanr + max_length;
2079  end1 = scan1 + BIGNUM_LENGTH(arg1);
2080  end2 = scan2 + BIGNUM_LENGTH(arg2);
2081
2082  carry2 = 1;
2083
2084  while (scanr < endr) {
2085    digit1 = (scan1 < end1) ? *scan1++ : 0;
2086    digit2 = (~((scan2 < end2) ? *scan2++ : 0) & BIGNUM_DIGIT_MASK)
2087             + carry2;
2088
2089    if (digit2 < BIGNUM_RADIX)
2090      carry2 = 0;
2091    else
2092      {
2093        digit2 = (digit2 - BIGNUM_RADIX);
2094        carry2 = 1;
2095      }
2096   
2097    *scanr++ = (op == bignum_and_op) ? digit1 & digit2 :
2098               (op == bignum_ior_op) ? digit1 | digit2 :
2099                                       digit1 ^ digit2;
2100  }
2101 
2102  if (neg_p)
2103    bignum_negate_magnitude(result);
2104
2105  return bignum_trim(result);
2106}
2107
2108static bignum_type
2109bignum_negneg_bitwise_op(int op, bignum_type arg1, bignum_type arg2)
2110{
2111  bignum_type result;
2112  bignum_length_type max_length;
2113
2114  bignum_digit_type *scan1, *end1, digit1, carry1;
2115  bignum_digit_type *scan2, *end2, digit2, carry2;
2116  bignum_digit_type *scanr, *endr;
2117
2118  char neg_p = op == bignum_and_op || op == bignum_ior_op;
2119
2120  max_length =  (BIGNUM_LENGTH(arg1) > BIGNUM_LENGTH(arg2))
2121               ? BIGNUM_LENGTH(arg1) + 1 : BIGNUM_LENGTH(arg2) + 1;
2122
2123  result = bignum_allocate(max_length, neg_p);
2124
2125  scanr = BIGNUM_START_PTR(result);
2126  scan1 = BIGNUM_START_PTR(arg1);
2127  scan2 = BIGNUM_START_PTR(arg2);
2128  endr = scanr + max_length;
2129  end1 = scan1 + BIGNUM_LENGTH(arg1);
2130  end2 = scan2 + BIGNUM_LENGTH(arg2);
2131
2132  carry1 = 1;
2133  carry2 = 1;
2134
2135  while (scanr < endr) {
2136    digit1 = (~((scan1 < end1) ? *scan1++ : 0) & BIGNUM_DIGIT_MASK) + carry1;
2137    digit2 = (~((scan2 < end2) ? *scan2++ : 0) & BIGNUM_DIGIT_MASK) + carry2;
2138
2139    if (digit1 < BIGNUM_RADIX)
2140      carry1 = 0;
2141    else
2142      {
2143        digit1 = (digit1 - BIGNUM_RADIX);
2144        carry1 = 1;
2145      }
2146   
2147    if (digit2 < BIGNUM_RADIX)
2148      carry2 = 0;
2149    else
2150      {
2151        digit2 = (digit2 - BIGNUM_RADIX);
2152        carry2 = 1;
2153      }
2154   
2155    *scanr++ = (op == bignum_and_op) ? digit1 & digit2 :
2156               (op == bignum_ior_op) ? digit1 | digit2 :
2157                                       digit1 ^ digit2;
2158  }
2159
2160  if (neg_p)
2161    bignum_negate_magnitude(result);
2162
2163  return bignum_trim(result);
2164}
2165
2166static void
2167bignum_negate_magnitude(bignum_type arg)
2168{
2169  bignum_digit_type *scan;
2170  bignum_digit_type *end;
2171  bignum_digit_type digit;
2172  bignum_digit_type carry;
2173
2174  scan = BIGNUM_START_PTR(arg);
2175  end = scan + BIGNUM_LENGTH(arg);
2176
2177  carry = 1;
2178
2179  while (scan < end) {
2180    digit = (~*scan & BIGNUM_DIGIT_MASK) + carry;
2181
2182    if (digit < BIGNUM_RADIX)
2183      carry = 0;
2184    else
2185      {
2186        digit = (digit - BIGNUM_RADIX);
2187        carry = 1;
2188      }
2189   
2190    *scan++ = digit;
2191  }
2192}
2193
2194#define BIGNUM_STR_BLOCK_SIZE 8
2195
2196/*
2197 * Contains s48_bignum_to_digit_stream.  We may want to separate that out
2198 * again in case of huge numbers; it may be more efficient to write
2199 * those straight to files without going through the string conversion.
2200 */
2201static void
2202big_to_string(C_word c, C_word self, C_word k, C_word value, C_word radix)
2203{
2204  char *buf, *index;
2205  C_word len, counter=0;
2206  char *characters = "0123456789abcdef";
2207  C_word *tmp, ret;
2208  bignum_type bignum;
2209
2210  bignum = big_of(value);
2211  len = BIGNUM_STR_BLOCK_SIZE;
2212  buf = C_malloc(len);
2213  if (buf == NULL) {
2214    fprintf(stderr, "out of memory - can not allocate string");
2215    exit(EXIT_FAILURE);
2216  }
2217  index = buf + len - 1;
2218
2219  if (BIGNUM_NEGATIVE_P(bignum)) {
2220    *index-- = '-';
2221    counter++;
2222  }
2223
2224  radix = C_unfix(radix);
2225 
2226  assert((radix > 1) && (radix <= BIGNUM_RADIX_ROOT));
2227  if (! (BIGNUM_ZERO_P (bignum)))
2228    {
2229      bignum_type working_copy = (bignum_copy (bignum));
2230      bignum_digit_type * start = (BIGNUM_START_PTR (working_copy));
2231      bignum_digit_type * scan = (start + (BIGNUM_LENGTH (working_copy)));
2232      bignum_digit_type digit;
2233      while (start < scan)
2234        {
2235          if ((scan[-1]) == 0) {
2236            scan -= 1;
2237          } else {
2238            digit = bignum_destructive_scale_down (working_copy, radix);
2239            *index-- = characters[digit];
2240            if (++counter == len) {
2241              char *newbuf;
2242              newbuf = C_malloc(len + BIGNUM_STR_BLOCK_SIZE);
2243              if (newbuf == NULL) {
2244                fprintf(stderr, "out of memory - can not allocate string");
2245                exit(EXIT_FAILURE);
2246              }
2247              C_memcpy(newbuf + BIGNUM_STR_BLOCK_SIZE, buf, len);
2248              C_free(buf);
2249              buf = newbuf;
2250              index = newbuf + BIGNUM_STR_BLOCK_SIZE - 1;
2251              len += BIGNUM_STR_BLOCK_SIZE;
2252            }
2253          }
2254        }
2255      BIGNUM_DEALLOCATE (working_copy);
2256    }
2257
2258  if (BIGNUM_NEGATIVE_P(bignum))
2259    *index-- = '-';
2260
2261  tmp = C_alloc(C_SIZEOF_STRING(counter));
2262  ret = C_string(&tmp, counter, index+1);
2263  C_free(buf);
2264  C_kontinue(k, ret);
2265}
2266
2267static void
2268digits_to_big(C_word c, C_word self, C_word k, C_word n,
2269              C_word start, C_word end, C_word radix, C_word negp)
2270{
2271  char *str;
2272  size_t n_digits;
2273  int negative_p = (negp != C_SCHEME_FALSE);
2274  int digit;
2275  int hash = 0;
2276
2277  str = C_c_string(n) + C_unfix(start);
2278  n_digits = C_unfix(end)-C_unfix(start);
2279  radix = C_unfix(radix);
2280  hash = /* abs(radix / 2) */ 0;
2281 
2282#define DIGIT_TO_INT(x)         \
2283  (((x) == '#') ? hash :        \
2284   (((x) >= (int)'a') ?((x) - (int)'a' + 10) : ((x) - (int)'0')))
2285
2286  assert((radix > 1) && (radix < BIGNUM_RADIX_ROOT));
2287  if (n_digits == 0)
2288    C_kontinue(k, C_SCHEME_FALSE);
2289  if (n_digits == 1)
2290    {
2291      digit = DIGIT_TO_INT(C_tolower((int)*str));
2292      if (digit >= radix || digit < 0)
2293        C_kontinue(k, C_SCHEME_FALSE);
2294      else
2295        C_return_bignum(k, bignum_digit_to_bignum(digit, negative_p));
2296    }
2297  {
2298    bignum_length_type length;
2299    {
2300      unsigned int radix_copy = radix;
2301      unsigned int log_radix = 0;
2302      while (radix_copy > 0)
2303        {
2304          radix_copy >>= 1;
2305          log_radix += 1;
2306        }
2307      /* This length will be at least as large as needed. */
2308      length = (BIGNUM_BITS_TO_DIGITS (n_digits * log_radix));
2309    }
2310    {
2311      bignum_type result = (bignum_allocate_zeroed (length, negative_p));
2312      while ((n_digits--) > 0)
2313        {
2314          digit = DIGIT_TO_INT(C_tolower((int)*str));
2315          str++;
2316          if (digit >= radix || digit < 0) {
2317            BIGNUM_DEALLOCATE(result);
2318            C_kontinue(k, C_SCHEME_FALSE);
2319          }
2320       
2321          bignum_destructive_scale_up (result, ((bignum_digit_type) radix));
2322          bignum_destructive_add (result, digit);
2323        }
2324      C_return_bignum (k, bignum_trim (result));
2325    }
2326  }
2327}
2328
2329/**
2330 * Below you will find the functions that have been refactored to
2331 * match the "core" style.
2332 *
2333 * Naming convention idea: Maybe name these fixnum ops differently,
2334 * _p_ for "promoting"?  OTOH, it may be safer and cleaner to rename
2335 * the old fixnum ops to have _fx_ to indicate they run in fixnum mode.
2336 */
2337static void bignum_negate_2(C_word c, C_word self, C_word new_big);
2338static void allocate_bignum_2(C_word c, C_word self, C_word bigvec);
2339
2340/* Eventually this will probably need to be integrated into C_2_plus. */
2341void C_ccall
2342C_u_2_fixnum_plus(C_word c, C_word self, C_word k, C_word x, C_word y)
2343{
2344  C_word z, ab[C_SIZEOF_BIGNUM(2)], *a = ab;
2345
2346  /* Exceptional situation: this will cause a real overflow */
2347  if (x == C_fix(C_MOST_NEGATIVE_FIXNUM) && y == C_fix(C_MOST_NEGATIVE_FIXNUM)) {
2348    C_kontinue(k, C_bignum2(&a, 1, 0, 2));
2349  } else {
2350    z = C_unfix(x) + C_unfix(y);
2351
2352    /* This code "knows" that both fixnums and bignums have 2 reserved bits */
2353    if(!C_fitsinfixnump(z)) {
2354      C_kontinue(k, C_bignum2(&a, (z < 0), labs(z) & (C_uword)BIGNUM_DIGIT_MASK, 1));
2355    } else {
2356      C_kontinue(k, C_fix(z));
2357    }
2358  }
2359}
2360
2361/* TODO: This should probably be renamed C_fixnum_negate to replace
2362 * what's in core.  Unfortunately, that one is allocating inline.
2363 * That one may be renamed to C_u_i_fixnum_negate() or some such.
2364 * TODO: Convert this to be an inline function and move to header?
2365 */
2366void C_ccall
2367C_u_fixnum_neg(C_word c, C_word self, C_word k, C_word x)
2368{
2369  /* Exceptional situation: this will cause an overflow to itself */
2370  if (x == C_fix(C_MOST_NEGATIVE_FIXNUM)) { /* C_fitsinfixnump(x) */
2371    C_word ab[C_SIZEOF_BIGNUM(2)], *a = ab;
2372    C_kontinue(k, C_bignum2(&a, 0, 0, 1));
2373  } else {
2374    C_kontinue(k, C_fix(-C_unfix(x)));
2375  }
2376}
2377
2378void C_ccall
2379C_fixnum_gcd(C_word c, C_word self, C_word k, C_word x, C_word y)
2380{
2381   C_word r;
2382   
2383   x = C_unfix(x);
2384   y = C_unfix(y);
2385   
2386   if (x < 0) x = -x;
2387   if (y < 0) y = -y;
2388   
2389   while(y != 0) {
2390     r = x % y;
2391     x = y;
2392     y = r;
2393   }
2394   C_kontinue(k, C_fix(x));
2395}
2396
2397void C_ccall
2398C_u_bignum_negate(C_word c, C_word self, C_word k, C_word x)
2399{
2400  C_word kab[C_SIZEOF_CLOSURE(3)], *ka = kab, k2, negp, size;
2401
2402  negp = C_i_not(C_u_i_bignum_negativep(x));
2403  k2 = C_closure(&ka, 3, (C_word)bignum_negate_2, k, x);
2404 
2405  size = C_u_i_bignum_size(x);
2406  C_allocate_bignum(3, (C_word)NULL, k2, size, negp, C_SCHEME_FALSE);
2407}
2408
2409static void
2410bignum_negate_2(C_word c, C_word self, C_word new_big)
2411{
2412  C_word k = C_block_item(self, 1),
2413         old_big = C_block_item(self, 2);
2414  C_memcpy(C_bignum_digits(new_big), C_bignum_digits(old_big),
2415           /* TODO: This is currently in bytes.  If we change the
2416            * representation that needs to change!
2417            * We subtract the size of the header, too.
2418            */
2419           C_header_size(C_internal_bignum(old_big))-C_wordstobytes(1));
2420  C_kontinue(k, new_big);
2421}
2422
2423void C_ccall
2424C_allocate_bignum(C_word c, C_word self, C_word k, C_word size, C_word negp, C_word initp)
2425{
2426  C_word kab[C_SIZEOF_CLOSURE(3)], *ka = kab, k2, init;
2427  k2 = C_closure(&ka, 3, (C_word)allocate_bignum_2, k, negp);
2428
2429  init = C_and(initp, C_make_character('\0'));
2430  C_allocate_vector(6, (C_word)NULL, k2,
2431                    C_bytes(C_fixnum_plus(size, C_fix(1))), /* Add header */
2432                    /* Byte vec, initialization, align at 8 bytes */
2433                    C_SCHEME_TRUE, init, C_SCHEME_FALSE);
2434}
2435
2436static void
2437allocate_bignum_2(C_word c, C_word self, C_word bigvec)
2438{
2439  C_word ab[C_SIZEOF_STRUCTURE(2)], *a = ab, bignum,
2440         k = C_block_item(self, 1),
2441         negp = C_truep(C_block_item(self, 2)),
2442         size = C_bytestowords(C_header_size(bigvec))-1;
2443
2444  C_word tagvec = CHICKEN_gc_root_ref(tags);
2445
2446  C_set_block_item(bigvec, 0, negp ? C_BIGNUM_HEADER_SIGN_BIT | size : size);
2447
2448  bignum = C_structure(&a, 2, C_block_item(tagvec, BIG_TAG), bigvec);
2449  C_kontinue(k, bignum);
2450}
Note: See TracBrowser for help on using the repository browser.