summaryrefslogtreecommitdiffstats
path: root/mpn_add.c
blob: d5c732a33577f7bc46c56d2ce251fa4df182aa7f (plain)
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
/* ----------------------------------------------------------------------- *
 *   
 *   Copyright 2007 H. Peter Anvin - All Rights Reserved
 *
 *   This program is free software; you can redistribute it and/or modify
 *   it under the terms of the GNU Lesser General Public License as
 *   published by the Free Software Foundation, Inc.,
 *   59 Temple Place Ste 330, Boston MA 02111-1307, USA; version 2.1,
 *   incorporated herein by reference.
 *
 * ----------------------------------------------------------------------- */

/*
 * mpn_add.c
 */

#include "mpnint.h"

struct mpn *mpn_addsub(struct mpn *s1, struct mpn *s2, int issub)
{
    int bits, len;
    struct mpn *d;
    int i;
    int s1minus = s1->minus;
    int s2minus = s2->minus ^ issub;
    int dminus;

    bits = s2->bits > s1->bits ? s2->bits+1 : s1->bits+1;
    len = (bits+MPN_ATOM_BITS-1)/MPN_ATOM_BITS;

    if (s1 == s2) {
	s1 = mpn_dupx(s1, len);
	mpn_free(s2);		/* Drop one reference */
    }

    if (s1minus == s2minus) {
	mpn_2atom_t c;
	/* Sign is the same, addition */

	if (s1->ref > s2->ref) {
	    struct mpn *t = s1;
	    s1 = s2;
	    s2 = t;
	}
	
	d = mpn_cow(s1, len);
	d->bits = bits;
	
	c = 0;
	for (i = 0; i < len; i++) {
	    c += d->num[i];
	    if (i < s2->len)
		c += s2->num[i];
	    
	    d->num[i] = c;
	    c >>= MPN_ATOM_BITS;
	}
	mpn_free(s2);

	/* d->minus already copied from s1->minus */
    } else {
	mpn_2satom_t c;		/* Relies on 2's complement with
				   proper right arithmetic shift... */
	/* We know the sign is different */
	dminus = 0;

	if (s2minus) {
	    struct mpn *t = s1;
	    s1 = s2;
	    s2 = t;
	    dminus = !dminus;
	}

	/* Now s1 is the addend and s2 is the subtrahend */
	if (mpn_cmp(s1, s2) < 0) {
	    /* addend is smaller than subtrahend, reverse and
	       invert minus */
	    struct mpn *t = s1;
	    s1 = s2;
	    s2 = t;
	    dminus = !dminus;
	}
	
	if (s1->ref > s2->ref) {
	    d = mpn_cow(s2, len);
	    
	    c = 0;
	    for (i = 0; i < len; i++) {
		c -= d->num[i];
		if (i < s1->len)
		    c += s1->num[i];
		
		d->num[i] = c;
		c >>= MPN_ATOM_BITS;
	    }
	    mpn_free(s1);
	} else {
	    d = mpn_cow(s1, len);
	    
	    c = 0;
	    for (i = 0; i < len; i++) {
		c += d->num[i];
		if (i < s2->len)
		    c -= s2->num[i];
		
		d->num[i] = c;
		c >>= MPN_ATOM_BITS;
	    }
	    mpn_free(s2);
	}

	d->minus = dminus;
    }

    _mpn_adjust_bits(d);

    return d;
}