1
2/*
3 *  M_APM  -  mapm_add.c
4 *
5 *  Copyright (C) 1999 - 2007   Michael C. Ring
6 *
7 *  Permission to use, copy, and distribute this software and its
8 *  documentation for any purpose with or without fee is hereby granted,
9 *  provided that the above copyright notice appear in all copies and
10 *  that both that copyright notice and this permission notice appear
11 *  in supporting documentation.
12 *
13 *  Permission to modify the software is granted. Permission to distribute
14 *  the modified code is granted. Modifications are to be distributed by
15 *  using the file 'license.txt' as a template to modify the file header.
16 *  'license.txt' is available in the official MAPM distribution.
17 *
18 *  This software is provided "as is" without express or implied warranty.
19 */
20
21/*
22 *      $Id: mapm_add.c,v 1.5 2007/12/03 01:33:39 mike Exp $
23 *
24 *      This file contains basic addition/subtraction functions
25 *
26 *      $Log: mapm_add.c,v $
27 *      Revision 1.5  2007/12/03 01:33:39  mike
28 *      Update license
29 *
30 *      Revision 1.4  2003/12/04 01:15:42  mike
31 *      redo math with 'borrow'
32 *
33 *      Revision 1.3  2002/11/03 22:03:31  mike
34 *      Updated function parameters to use the modern style
35 *
36 *      Revision 1.2  2001/07/16 18:59:25  mike
37 *      add function M_free_all_add
38 *
39 *      Revision 1.1  1999/05/10 20:56:31  mike
40 *      Initial revision
41 */
42
43#include "m_apm_lc.h"
44
45static	M_APM	M_work1 = NULL;
46static	M_APM	M_work2 = NULL;
47static	int	M_add_firsttime = TRUE;
48
49/****************************************************************************/
50void	M_free_all_add()
51{
52if (M_add_firsttime == FALSE)
53  {
54   m_apm_free(M_work1);
55   m_apm_free(M_work2);
56   M_add_firsttime = TRUE;
57  }
58}
59/****************************************************************************/
60void	m_apm_add(M_APM r, M_APM a, M_APM b)
61{
62int	j, carry, sign, aexp, bexp, adigits, bdigits;
63
64if (M_add_firsttime)
65  {
66   M_add_firsttime = FALSE;
67   M_work1 = m_apm_init();
68   M_work2 = m_apm_init();
69  }
70
71if (a->m_apm_sign == 0)
72  {
73   m_apm_copy(r,b);
74   return;
75  }
76
77if (b->m_apm_sign == 0)
78  {
79   m_apm_copy(r,a);
80   return;
81  }
82
83if (a->m_apm_sign == 1 && b->m_apm_sign == -1)
84  {
85   b->m_apm_sign = 1;
86   m_apm_subtract(r,a,b);
87   b->m_apm_sign = -1;
88   return;
89  }
90
91if (a->m_apm_sign == -1 && b->m_apm_sign == 1)
92  {
93   a->m_apm_sign = 1;
94   m_apm_subtract(r,b,a);
95   a->m_apm_sign = -1;
96   return;
97  }
98
99sign = a->m_apm_sign;         /* signs are the same, result will be same */
100
101aexp = a->m_apm_exponent;
102bexp = b->m_apm_exponent;
103
104m_apm_copy(M_work1, a);
105m_apm_copy(M_work2, b);
106
107/*
108 *  scale by at least 1 factor of 10 in case the MSB carrys
109 */
110
111if (aexp == bexp)
112  {
113   M_apm_scale(M_work1, 2);   /* shift 2 digits == 1 byte for efficiency */
114   M_apm_scale(M_work2, 2);
115  }
116else
117  {
118   if (aexp > bexp)
119     {
120      M_apm_scale(M_work1, 2);
121      M_apm_scale(M_work2, (aexp + 2 - bexp));
122     }
123   else            /*  aexp < bexp  */
124     {
125      M_apm_scale(M_work2, 2);
126      M_apm_scale(M_work1, (bexp + 2 - aexp));
127     }
128  }
129
130adigits = M_work1->m_apm_datalength;
131bdigits = M_work2->m_apm_datalength;
132
133if (adigits >= bdigits)
134  {
135   m_apm_copy(r, M_work1);
136   j = (bdigits + 1) >> 1;
137   carry = 0;
138
139   while (TRUE)
140     {
141      j--;
142      r->m_apm_data[j] += carry + M_work2->m_apm_data[j];
143
144      if (r->m_apm_data[j] >= 100)
145        {
146         r->m_apm_data[j] -= 100;
147	 carry = 1;
148	}
149      else
150         carry = 0;
151
152      if (j == 0)
153        break;
154     }
155  }
156else
157  {
158   m_apm_copy(r, M_work2);
159   j = (adigits + 1) >> 1;
160   carry = 0;
161
162   while (TRUE)
163     {
164      j--;
165      r->m_apm_data[j] += carry + M_work1->m_apm_data[j];
166
167      if (r->m_apm_data[j] >= 100)
168        {
169         r->m_apm_data[j] -= 100;
170	 carry = 1;
171	}
172      else
173         carry = 0;
174
175      if (j == 0)
176        break;
177     }
178  }
179
180r->m_apm_sign = sign;
181
182M_apm_normalize(r);
183}
184/****************************************************************************/
185void	m_apm_subtract(M_APM r, M_APM a, M_APM b)
186{
187int	itmp, j, flag, icompare, sign, aexp, bexp,
188	borrow, adigits, bdigits;
189
190if (M_add_firsttime)
191  {
192   M_add_firsttime = FALSE;
193   M_work1 = m_apm_init();
194   M_work2 = m_apm_init();
195  }
196
197if (b->m_apm_sign == 0)
198  {
199   m_apm_copy(r,a);
200   return;
201  }
202
203if (a->m_apm_sign == 0)
204  {
205   m_apm_copy(r,b);
206   r->m_apm_sign = -(r->m_apm_sign);
207   return;
208  }
209
210if (a->m_apm_sign == 1 && b->m_apm_sign == -1)
211  {
212   b->m_apm_sign = 1;
213   m_apm_add(r,a,b);
214   b->m_apm_sign = -1;
215   return;
216  }
217
218if (a->m_apm_sign == -1 && b->m_apm_sign == 1)
219  {
220   b->m_apm_sign = -1;
221   m_apm_add(r,a,b);
222   b->m_apm_sign = 1;
223   return;
224  }
225
226/* now, the signs are the same  */
227/* make a positive working copy */
228
229m_apm_absolute_value(M_work1, a);
230m_apm_absolute_value(M_work2, b);
231
232/* are they the same??  if so, the result is zero */
233
234if ((icompare = m_apm_compare(M_work1, M_work2)) == 0)
235  {
236   M_set_to_zero(r);
237   return;
238  }
239
240if (icompare == 1)             /*  |a| > |b|  (do A-B)  */
241  {
242   flag = TRUE;
243   sign = a->m_apm_sign;
244  }
245else                           /*  |b| > |a|  (do B-A)  */
246  {
247   flag = FALSE;
248   sign = -(a->m_apm_sign);
249  }
250
251aexp = M_work1->m_apm_exponent;
252bexp = M_work2->m_apm_exponent;
253
254if (aexp > bexp)
255  M_apm_scale(M_work2, (aexp - bexp));
256
257if (aexp < bexp)
258  M_apm_scale(M_work1, (bexp - aexp));
259
260adigits = M_work1->m_apm_datalength;
261bdigits = M_work2->m_apm_datalength;
262
263if (adigits > bdigits)
264  M_apm_pad(M_work2, adigits);
265
266if (adigits < bdigits)
267  M_apm_pad(M_work1, bdigits);
268
269if (flag)		/* perform A-B,  M_work1 - M_work2 */
270  {
271   m_apm_copy(r, M_work1);
272   j = (r->m_apm_datalength + 1) >> 1;
273   borrow = 0;
274
275   while (TRUE)
276     {
277      j--;
278      itmp = (int)r->m_apm_data[j] - ((int)M_work2->m_apm_data[j] + borrow);
279
280      if (itmp >= 0)
281        {
282         r->m_apm_data[j] = (UCHAR)itmp;
283	 borrow = 0;
284        }
285      else
286        {
287         r->m_apm_data[j] = (UCHAR)(100 + itmp);
288	 borrow = 1;
289	}
290
291      if (j == 0)
292        break;
293     }
294  }
295else   		/* perform B-A,  M_work2 - M_work1 */
296  {
297   m_apm_copy(r, M_work2);
298   j = (r->m_apm_datalength + 1) >> 1;
299   borrow = 0;
300
301   while (TRUE)
302     {
303      j--;
304      itmp = (int)r->m_apm_data[j] - ((int)M_work1->m_apm_data[j] + borrow);
305
306      if (itmp >= 0)
307        {
308         r->m_apm_data[j] = (UCHAR)itmp;
309	 borrow = 0;
310        }
311      else
312        {
313         r->m_apm_data[j] = (UCHAR)(100 + itmp);
314	 borrow = 1;
315	}
316
317      if (j == 0)
318        break;
319     }
320  }
321
322r->m_apm_sign = sign;
323
324M_apm_normalize(r);
325}
326/****************************************************************************/
327