ap.c
资源名称:c.rar [点击查看]
上传用户:shmaik
上传日期:2014-06-01
资源大小:45093k
文件大小:9k
源码类别:

VC书籍

开发平台:

C/C++

  1. static char rcsid[] = "$Id: H:/drh/idioms/book/RCS/ap.doc,v 1.11 1996/06/26 23:02:01 drh Exp $";
  2. #include <ctype.h>
  3. #include <limits.h>
  4. #include <stdlib.h>
  5. #include <string.h>
  6. #include "assert.h"
  7. #include "ap.h"
  8. #include "fmt.h"
  9. #include "xp.h"
  10. #include "mem.h"
  11. #define T AP_T
  12. struct T {
  13. int sign;
  14. int ndigits;
  15. int size;
  16. XP_T digits;
  17. };
  18. #define iszero(x) ((x)->ndigits==1 && (x)->digits[0]==0)
  19. #define maxdigits(x,y) ((x)->ndigits > (y)->ndigits ? 
  20. (x)->ndigits : (y)->ndigits)
  21. #define isone(x) ((x)->ndigits==1 && (x)->digits[0]==1)
  22. static T normalize(T z, int n);
  23. static int cmp(T x, T y);
  24. static T mk(int size) {
  25. T z = CALLOC(1, sizeof (*z) + size);
  26. assert(size > 0);
  27. z->sign = 1;
  28. z->size = size;
  29. z->ndigits = 1;
  30. z->digits = (XP_T)(z + 1);
  31. return z;
  32. }
  33. static T set(T z, long int n) {
  34. if (n == LONG_MIN)
  35. XP_fromint(z->size, z->digits, LONG_MAX + 1UL);
  36. else if (n < 0)
  37. XP_fromint(z->size, z->digits, -n);
  38. else
  39. XP_fromint(z->size, z->digits, n);
  40. z->sign = n < 0 ? -1 : 1;
  41. return normalize(z, z->size);
  42. }
  43. static T normalize(T z, int n) {
  44. z->ndigits = XP_length(n, z->digits);
  45. return z;
  46. }
  47. static T add(T z, T x, T y) {
  48. int n = y->ndigits;
  49. if (x->ndigits < n)
  50. return add(z, y, x);
  51. else if (x->ndigits > n) {
  52. int carry = XP_add(n, z->digits, x->digits,
  53. y->digits, 0);
  54. z->digits[z->size-1] = XP_sum(x->ndigits - n,
  55. &z->digits[n], &x->digits[n], carry);
  56. } else
  57. z->digits[n] = XP_add(n, z->digits, x->digits,
  58. y->digits, 0);
  59. return normalize(z, z->size);
  60. }
  61. static T sub(T z, T x, T y) {
  62. int borrow, n = y->ndigits;
  63. borrow = XP_sub(n, z->digits, x->digits,
  64. y->digits, 0);
  65. if (x->ndigits > n)
  66. borrow = XP_diff(x->ndigits - n, &z->digits[n],
  67. &x->digits[n], borrow);
  68. assert(borrow == 0);
  69. return normalize(z, z->size);
  70. }
  71. static T mulmod(T x, T y, T p) {
  72. T z, xy = AP_mul(x, y);
  73. z = AP_mod(xy, p);
  74. AP_free(&xy);
  75. return z;
  76. }
  77. static int cmp(T x, T y) {
  78. if (x->ndigits != y->ndigits)
  79. return x->ndigits - y->ndigits;
  80. else
  81. return XP_cmp(x->ndigits, x->digits, y->digits);
  82. }
  83. T AP_new(long int n) {
  84. return set(mk(sizeof (long int)), n);
  85. }
  86. void AP_free(T *z) {
  87. assert(z && *z);
  88. FREE(*z);
  89. }
  90. T AP_neg(T x) {
  91. T z;
  92. assert(x);
  93. z = mk(x->ndigits);
  94. memcpy(z->digits, x->digits, x->ndigits);
  95. z->ndigits = x->ndigits;
  96. z->sign = iszero(z) ? 1 : -x->sign;
  97. return z;
  98. }
  99. T AP_mul(T x, T y) {
  100. T z;
  101. assert(x);
  102. assert(y);
  103. z = mk(x->ndigits + y->ndigits);
  104. XP_mul(z->digits, x->ndigits, x->digits, y->ndigits,
  105. y->digits);
  106. normalize(z, z->size);
  107. z->sign = iszero(z)
  108. || ((x->sign^y->sign) == 0) ? 1 : -1;
  109. return z;
  110. }
  111. T AP_add(T x, T y) {
  112. T z;
  113. assert(x);
  114. assert(y);
  115. if (((x->sign^y->sign) == 0)) {
  116. z = add(mk(maxdigits(x,y) + 1), x, y);
  117. z->sign = iszero(z) ? 1 : x->sign;
  118. } else
  119. if (cmp(x, y) > 0) {
  120. z = sub(mk(x->ndigits), x, y);
  121. z->sign = iszero(z) ? 1 : x->sign;
  122. }
  123. else {
  124. z = sub(mk(y->ndigits), y, x);
  125. z->sign = iszero(z) ? 1 : -x->sign;
  126. }
  127. return z;
  128. }
  129. T AP_sub(T x, T y) {
  130. T z;
  131. assert(x);
  132. assert(y);
  133. if (!((x->sign^y->sign) == 0)) {
  134. z = add(mk(maxdigits(x,y) + 1), x, y);
  135. z->sign = iszero(z) ? 1 : x->sign;
  136. } else
  137. if (cmp(x, y) > 0) {
  138. z = sub(mk(x->ndigits), x, y);
  139. z->sign = iszero(z) ? 1 : x->sign;
  140. } else {
  141. z = sub(mk(y->ndigits), y, x);
  142. z->sign = iszero(z) ? 1 : -x->sign;
  143. }
  144. return z;
  145. }
  146. T AP_div(T x, T y) {
  147. T q, r;
  148. assert(x);
  149. assert(y);
  150. assert(!iszero(y));
  151. q = mk(x->ndigits);
  152. r = mk(y->ndigits);
  153. {
  154. XP_T tmp = ALLOC(x->ndigits + y->ndigits + 2);
  155. XP_div(x->ndigits, q->digits, x->digits,
  156. y->ndigits, y->digits, r->digits, tmp);
  157. FREE(tmp);
  158. }
  159. normalize(q, q->size);
  160. normalize(r, r->size);
  161. q->sign = iszero(q)
  162. || ((x->sign^y->sign) == 0) ? 1 : -1;
  163. if (!((x->sign^y->sign) == 0) && !iszero(r)) {
  164. int carry = XP_sum(q->size, q->digits,
  165. q->digits, 1);
  166. assert(carry == 0);
  167. normalize(q, q->size);
  168. }
  169. AP_free(&r);
  170. return q;
  171. }
  172. T AP_mod(T x, T y) {
  173. T q, r;
  174. assert(x);
  175. assert(y);
  176. assert(!iszero(y));
  177. q = mk(x->ndigits);
  178. r = mk(y->ndigits);
  179. {
  180. XP_T tmp = ALLOC(x->ndigits + y->ndigits + 2);
  181. XP_div(x->ndigits, q->digits, x->digits,
  182. y->ndigits, y->digits, r->digits, tmp);
  183. FREE(tmp);
  184. }
  185. normalize(q, q->size);
  186. normalize(r, r->size);
  187. q->sign = iszero(q)
  188. || ((x->sign^y->sign) == 0) ? 1 : -1;
  189. if (!((x->sign^y->sign) == 0) && !iszero(r)) {
  190. int borrow = XP_sub(r->size, r->digits,
  191. y->digits, r->digits, 0);
  192. assert(borrow == 0);
  193. normalize(r, r->size);
  194. }
  195. AP_free(&q);
  196. return r;
  197. }
  198. T AP_pow(T x, T y, T p) {
  199. T z;
  200. assert(x);
  201. assert(y);
  202. assert(y->sign == 1);
  203. assert(!p || p->sign==1 && !iszero(p) && !isone(p));
  204. if (iszero(x))
  205. return AP_new(0);
  206. if (iszero(y))
  207. return AP_new(1);
  208. if (isone(x))
  209. return AP_new((((y)->digits[0]&1) == 0) ? 1 : x->sign);
  210. if (p)
  211. if (isone(y))
  212. z = AP_mod(x, p);
  213. else {
  214. T y2 = AP_rshift(y, 1), t = AP_pow(x, y2, p);
  215. z = mulmod(t, t, p);
  216. AP_free(&y2);
  217. AP_free(&t);
  218. if (!(((y)->digits[0]&1) == 0)) {
  219. z = mulmod(y2 = AP_mod(x, p), t = z, p);
  220. AP_free(&y2);
  221. AP_free(&t);
  222. }
  223. }
  224. else
  225. if (isone(y))
  226. z = AP_addi(x, 0);
  227. else {
  228. T y2 = AP_rshift(y, 1), t = AP_pow(x, y2, NULL);
  229. z = AP_mul(t, t);
  230. AP_free(&y2);
  231. AP_free(&t);
  232. if (!(((y)->digits[0]&1) == 0)) {
  233. z = AP_mul(x, t = z);
  234. AP_free(&t);
  235. }
  236. }
  237. return z;
  238. }
  239. int AP_cmp(T x, T y) {
  240. assert(x);
  241. assert(y);
  242. if (!((x->sign^y->sign) == 0))
  243. return x->sign;
  244. else if (x->sign == 1)
  245. return cmp(x, y);
  246. else
  247. return cmp(y, x);
  248. }
  249. T AP_addi(T x, long int y) {
  250. unsigned char d[sizeof (unsigned long)];
  251. struct T t;
  252. t.size = sizeof d;
  253. t.digits = d;
  254. return AP_add(x, set(&t, y));
  255. }
  256. T AP_subi(T x, long int y) {
  257. unsigned char d[sizeof (unsigned long)];
  258. struct T t;
  259. t.size = sizeof d;
  260. t.digits = d;
  261. return AP_sub(x, set(&t, y));
  262. }
  263. T AP_muli(T x, long int y) {
  264. unsigned char d[sizeof (unsigned long)];
  265. struct T t;
  266. t.size = sizeof d;
  267. t.digits = d;
  268. return AP_mul(x, set(&t, y));
  269. }
  270. T AP_divi(T x, long int y) {
  271. unsigned char d[sizeof (unsigned long)];
  272. struct T t;
  273. t.size = sizeof d;
  274. t.digits = d;
  275. return AP_div(x, set(&t, y));
  276. }
  277. int AP_cmpi(T x, long int y) {
  278. unsigned char d[sizeof (unsigned long)];
  279. struct T t;
  280. t.size = sizeof d;
  281. t.digits = d;
  282. return AP_cmp(x, set(&t, y));
  283. }
  284. long int AP_modi(T x, long int y) {
  285. long int rem;
  286. T r;
  287. unsigned char d[sizeof (unsigned long)];
  288. struct T t;
  289. t.size = sizeof d;
  290. t.digits = d;
  291. r = AP_mod(x, set(&t, y));
  292. rem = XP_toint(r->ndigits, r->digits);
  293. AP_free(&r);
  294. return rem;
  295. }
  296. T AP_lshift(T x, int s) {
  297. T z;
  298. assert(x);
  299. assert(s >= 0);
  300. z = mk(x->ndigits + ((s+7)&~7)/8);
  301. XP_lshift(z->size, z->digits, x->ndigits,
  302. x->digits, s, 0);
  303. z->sign = x->sign;
  304. return normalize(z, z->size);
  305. }
  306. T AP_rshift(T x, int s) {
  307. assert(x);
  308. assert(s >= 0);
  309. if (s >= 8*x->ndigits)
  310. return AP_new(0);
  311. else {
  312. T z = mk(x->ndigits - s/8);
  313. XP_rshift(z->size, z->digits, x->ndigits,
  314. x->digits, s, 0);
  315. normalize(z, z->size);
  316. z->sign = iszero(z) ? 1 : x->sign;
  317. return z;
  318. }
  319. }
  320. long int AP_toint(T x) {
  321. unsigned long u;
  322. assert(x);
  323. u = XP_toint(x->ndigits, x->digits)%(LONG_MAX + 1UL);
  324. if (x->sign == -1)
  325. return -(long)u;
  326. else
  327. return  (long)u;
  328. }
  329. T AP_fromstr(const char *str, int base, char **end) {
  330. T z;
  331. const char *p = str;
  332. char *endp, sign = '';
  333. int carry;
  334. assert(p);
  335. assert(base >= 2 && base <= 36);
  336. while (*p && isspace(*p))
  337. p++;
  338. if (*p == '-' || *p == '+')
  339. sign = *p++;
  340. {
  341. const char *start;
  342. int k, n = 0;
  343. for ( ; *p == '0' && p[1] == '0'; p++)
  344. ;
  345. start = p;
  346. for ( ; (  '0' <= *p && *p <= '9' && *p < '0' + base
  347. || 'a' <= *p && *p <= 'z' && *p < 'a' + base - 10
  348. || 'A' <= *p && *p <= 'Z' && *p < 'A' + base - 10); p++)
  349. n++;
  350. for (k = 1; (1<<k) < base; k++)
  351. ;
  352. z = mk(((k*n + 7)&~7)/8);
  353. p = start;
  354. }
  355. carry = XP_fromstr(z->size, z->digits, p,
  356. base, &endp);
  357. assert(carry == 0);
  358. normalize(z, z->size);
  359. if (endp == p) {
  360. endp = (char *)str;
  361. z = AP_new(0);
  362. } else
  363. z->sign = iszero(z) || sign != '-' ? 1 : -1;
  364. if (end)
  365. *end = (char *)endp;
  366. return z;
  367. }
  368. char *AP_tostr(char *str, int size, int base, T x) {
  369. XP_T q;
  370. assert(x);
  371. assert(base >= 2 && base <= 36);
  372. assert(str == NULL || size > 1);
  373. if (str == NULL) {
  374. {
  375. int k;
  376. for (k = 5; (1<<k) > base; k--)
  377. ;
  378. size = (8*x->ndigits)/k + 1 + 1;
  379. if (x->sign == 1)
  380. size++;
  381. }
  382. str = ALLOC(size);
  383. }
  384. q = ALLOC(x->ndigits);
  385. memcpy(q, x->digits, x->ndigits);
  386. if (x->sign == -1) {
  387. str[0] = '-';
  388. XP_tostr(str + 1, size - 1, base, x->ndigits, q);
  389. } else
  390. XP_tostr(str, size, base, x->ndigits, q);
  391. FREE(q);
  392. return str;
  393. }
  394. void AP_fmt(int code, va_list *app,
  395. int put(int c, void *cl), void *cl,
  396. unsigned char flags[], int width, int precision) {
  397. T x;
  398. char *buf;
  399. assert(app && flags);
  400. x = va_arg(*app, T);
  401. assert(x);
  402. buf = AP_tostr(NULL, 0, 10, x);
  403. Fmt_putd(buf, strlen(buf), put, cl, flags,
  404. width, precision);
  405. FREE(buf);
  406. }