mp 10 KB

123456789101112131415161718192021222324252627282930313233343536373839404142434445464748495051525354555657585960616263646566676869707172737475767778798081828384858687888990919293949596979899100101102103104105106107108109110111112113114115116117118119120121122123124125126127128129130131132133134135136137138139140141142143144145146147148149150151152153154155156157158159160161162163164165166167168169170171172173174175176177178179180181182183184185186187188189190191192193194195196197198199200201202203204205206207208209210211212213214215216217218219220221222223224225226227228229230231232233234235236237238239240241242243244245246247248249250251252253254255256257258259260261262263264265266267268269270271272273274275276277278279280281282283284285286287288289290291292293294295296297298299300301302303304305306307308309310311312313314315316317318319320321322323324325326327328329330331332333334335336337338339340341342343344345346347348349350351352353354355356357358359360361362363364365366367368369370371372373374375376377378379380381382383384385386387388389390391392393394395396397398399400401402403404405406407408409410411412413414415416417418419420421422423424425426427428429430431432433434435436437438439440441442443444445446447448449450451452453454455456457458459460461462463464465466467468469470471472473474475476477478479480481482483484485486487488489490491492493494495496497498499500501502503504505506507508509510511512513514515516517518519520521522523524525526527528529530531532533534535536537538539540541542543544545546547548549550551552553554555556557558559560561562563564565566567568569570571572573574575576577
  1. .TH MP 2
  2. .SH NAME
  3. mpsetminbits, mpnew, mpfree, mpbits, mpnorm, mpcopy, mpassign, mprand, strtomp, mpfmt,mptoa, betomp, mptobe, letomp, mptole, mptoui, uitomp, mptoi, itomp, uvtomp, mptouv, vtomp, mptov, mpdigdiv, mpadd, mpsub, mpleft, mpright, mpmul, mpexp, mpmod, mpdiv, mpcmp, mpextendedgcd, mpinvert, mpsignif, mplowbits0, mpvecdigmuladd, mpvecdigmulsub, mpvecadd, mpvecsub, mpveccmp, mpvecmul, mpmagcmp, mpmagadd, mpmagsub, crtpre, crtin, crtout, crtprefree, crtresfree \- extended precision arithmetic
  4. .SH SYNOPSIS
  5. .B #include <u.h>
  6. .br
  7. .B #include <libc.h>
  8. .br
  9. .B #include <mp.h>
  10. .PP
  11. .B
  12. mpint* mpnew(int n)
  13. .PP
  14. .B
  15. void mpfree(mpint *b)
  16. .PP
  17. .B
  18. void mpsetminbits(int n)
  19. .PP
  20. .B
  21. void mpbits(mpint *b, int n)
  22. .PP
  23. .B
  24. void mpnorm(mpint *b)
  25. .PP
  26. .B
  27. mpint* mpcopy(mpint *b)
  28. .PP
  29. .B
  30. void mpassign(mpint *old, mpint *new)
  31. .PP
  32. .B
  33. mpint* mprand(int bits, void (*gen)(uchar*, int), mpint *b)
  34. .PP
  35. .B
  36. mpint* strtomp(char *buf, char **rptr, int base, mpint *b)
  37. .PP
  38. .B
  39. char* mptoa(mpint *b, int base, char *buf, int blen)
  40. .PP
  41. .B
  42. int mpfmt(Fmt*)
  43. .PP
  44. .B
  45. mpint* betomp(uchar *buf, uint blen, mpint *b)
  46. .PP
  47. .B
  48. int mptobe(mpint *b, uchar *buf, uint blen, uchar **bufp)
  49. .PP
  50. .B
  51. mpint* letomp(uchar *buf, uint blen, mpint *b)
  52. .PP
  53. .B
  54. int mptole(mpint *b, uchar *buf, uint blen, uchar **bufp)
  55. .PP
  56. .B
  57. uint mptoui(mpint*)
  58. .PP
  59. .B
  60. mpint* uitomp(uint, mpint*)
  61. .PP
  62. .B
  63. int mptoi(mpint*)
  64. .PP
  65. .B
  66. mpint* itomp(int, mpint*)
  67. .PP
  68. .B
  69. mpint* vtomp(vlong, mpint*)
  70. .PP
  71. .B
  72. vlong mptov(mpint*)
  73. .PP
  74. .B
  75. mpint* uvtomp(uvlong, mpint*)
  76. .PP
  77. .B
  78. uvlong mptouv(mpint*)
  79. .PP
  80. .B
  81. void mpadd(mpint *b1, mpint *b2, mpint *sum)
  82. .PP
  83. .B
  84. void mpmagadd(mpint *b1, mpint *b2, mpint *sum)
  85. .PP
  86. .B
  87. void mpsub(mpint *b1, mpint *b2, mpint *diff)
  88. .PP
  89. .B
  90. void mpmagsub(mpint *b1, mpint *b2, mpint *diff)
  91. .PP
  92. .B
  93. void mpleft(mpint *b, int shift, mpint *res)
  94. .PP
  95. .B
  96. void mpright(mpint *b, int shift, mpint *res)
  97. .PP
  98. .B
  99. void mpmul(mpint *b1, mpint *b2, mpint *prod)
  100. .PP
  101. .B
  102. void mpexp(mpint *b, mpint *e, mpint *m, mpint *res)
  103. .PP
  104. .B
  105. void mpmod(mpint *b, mpint *m, mpint *remainder)
  106. .PP
  107. .B
  108. void mpdiv(mpint *dividend, mpint *divisor, mpint *quotient, mpint *remainder)
  109. .PP
  110. .B
  111. int mpcmp(mpint *b1, mpint *b2)
  112. .PP
  113. .B
  114. int mpmagcmp(mpint *b1, mpint *b2)
  115. .PP
  116. .B
  117. void mpextendedgcd(mpint *a, mpint *b, mpint *d, mpint *x, mpint *y)
  118. .PP
  119. .B
  120. void mpinvert(mpint *b, mpint *m, mpint *res)
  121. .PP
  122. .B
  123. int mpsignif(mpint *b)
  124. .PP
  125. .B
  126. int mplowbits0(mpint *b)
  127. .PP
  128. .B
  129. void mpdigdiv(mpdigit *dividend, mpdigit divisor, mpdigit *quotient)
  130. .PP
  131. .B
  132. void mpvecadd(mpdigit *a, int alen, mpdigit *b, int blen, mpdigit *sum)
  133. .PP
  134. .B
  135. void mpvecsub(mpdigit *a, int alen, mpdigit *b, int blen, mpdigit *diff)
  136. .PP
  137. .B
  138. void mpvecdigmuladd(mpdigit *b, int n, mpdigit m, mpdigit *p)
  139. .PP
  140. .B
  141. int mpvecdigmulsub(mpdigit *b, int n, mpdigit m, mpdigit *p)
  142. .PP
  143. .B
  144. void mpvecmul(mpdigit *a, int alen, mpdigit *b, int blen, mpdigit *p)
  145. .PP
  146. .B
  147. int mpveccmp(mpdigit *a, int alen, mpdigit *b, int blen)
  148. .PP
  149. .B
  150. CRTpre* crtpre(int nfactors, mpint **factors)
  151. .PP
  152. .B
  153. CRTres* crtin(CRTpre *crt, mpint *x)
  154. .PP
  155. .B
  156. void crtout(CRTpre *crt, CRTres *r, mpint *x)
  157. .PP
  158. .B
  159. void crtprefree(CRTpre *cre)
  160. .PP
  161. .B
  162. void crtresfree(CRTres *res)
  163. .PP
  164. .B
  165. mpint *mpzero, *mpone, *mptwo
  166. .SH DESCRIPTION
  167. .PP
  168. These routines perform extended precision integer arithmetic.
  169. The basic type is
  170. .BR mpint ,
  171. which points to an array of
  172. .BR mpdigit s,
  173. stored in little-endian order:
  174. .sp
  175. .EX
  176. typedef struct mpint mpint;
  177. struct mpint
  178. {
  179. int sign; /* +1 or -1 */
  180. int size; /* allocated digits */
  181. int top; /* significant digits */
  182. mpdigit *p;
  183. char flags;
  184. };
  185. .EE
  186. .PP
  187. The sign of 0 is +1.
  188. .PP
  189. The size of
  190. .B mpdigit
  191. is architecture-dependent and defined in
  192. .BR /$cputype/include/u.h .
  193. .BR Mpint s
  194. are dynamically allocated and must be explicitly freed. Operations
  195. grow the array of digits as needed.
  196. .PP
  197. In general, the result parameters are last in the
  198. argument list.
  199. .PP
  200. Routines that return an
  201. .B mpint
  202. will allocate the
  203. .B mpint
  204. if the result parameter is
  205. .BR nil .
  206. This includes
  207. .IR strtomp ,
  208. .IR itomp ,
  209. .IR uitomp ,
  210. and
  211. .IR btomp .
  212. These functions, in addition to
  213. .I mpnew
  214. and
  215. .IR mpcopy ,
  216. will return
  217. .B nil
  218. if the allocation fails.
  219. .PP
  220. Input and result parameters may point to the same
  221. .BR mpint .
  222. The routines check and copy where necessary.
  223. .PP
  224. .I Mpnew
  225. creates an
  226. .B mpint
  227. with an initial allocation of
  228. .I n
  229. bits.
  230. If
  231. .I n
  232. is zero, the allocation will be whatever was specified in the
  233. last call to
  234. .I mpsetminbits
  235. or to the initial value, 1056.
  236. .I Mpfree
  237. frees an
  238. .BR mpint .
  239. .I Mpbits
  240. grows the allocation of
  241. .I b
  242. to fit at least
  243. .I n
  244. bits. If
  245. .B b->top
  246. doesn't cover
  247. .I n
  248. bits it increases it to do so.
  249. Unless you are writing new basic operations, you
  250. can restrict yourself to
  251. .B mpnew(0)
  252. and
  253. .BR mpfree(b) .
  254. .PP
  255. .I Mpnorm
  256. normalizes the representation by trimming any high order zero
  257. digits. All routines except
  258. .B mpbits
  259. return normalized results.
  260. .PP
  261. .I Mpcopy
  262. creates a new
  263. .B mpint
  264. with the same value as
  265. .I b
  266. while
  267. .I mpassign
  268. sets the value of
  269. .I new
  270. to be that of
  271. .IR old .
  272. .PP
  273. .I Mprand
  274. creates an
  275. .I n
  276. bit random number using the generator
  277. .IR gen .
  278. .I Gen
  279. takes a pointer to a string of uchar's and the number
  280. to fill in.
  281. .PP
  282. .I Strtomp
  283. and
  284. .I mptoa
  285. convert between
  286. .SM ASCII
  287. and
  288. .B mpint
  289. representations using the base indicated.
  290. Only the bases 10, 16, 32, and 64 are
  291. supported. Anything else defaults to 16.
  292. .IR Strtomp
  293. skips any leading spaces or tabs.
  294. .IR Strtomp 's
  295. scan stops when encountering a digit not valid in the
  296. base. If
  297. .I rptr
  298. is not zero,
  299. .I *rptr
  300. is set to point to the character immediately after the
  301. string converted.
  302. If the parse pterminates before any digits are found,
  303. .I strtomp
  304. return
  305. .BR nil .
  306. .I Mptoa
  307. returns a pointer to the filled buffer.
  308. If the parameter
  309. .I buf
  310. is
  311. .BR nil ,
  312. the buffer is allocated.
  313. .I Mpfmt
  314. can be used with
  315. .IR fmtinstall (2)
  316. and
  317. .IR print (2)
  318. to print hexadecimal representations of
  319. .BR mpint s.
  320. .PP
  321. .I Mptobe
  322. and
  323. .I mptole
  324. convert an
  325. .I mpint
  326. to a byte array. The former creates a big endian representation,
  327. the latter a little endian one.
  328. If the destination
  329. .I buf
  330. is not
  331. .BR nil ,
  332. it specifies the buffer of length
  333. .I blen
  334. for the result. If the representation
  335. is less than
  336. .I blen
  337. bytes, the rest of the buffer is zero filled.
  338. If
  339. .I buf
  340. is
  341. .BR nil ,
  342. then a buffer is allocated and a pointer to it is
  343. deposited in the location pointed to by
  344. .IR bufp .
  345. Sign is ignored in these conversions, i.e., the byte
  346. array version is always positive.
  347. .PP
  348. .IR Betomp ,
  349. and
  350. .I letomp
  351. convert from a big or little endian byte array at
  352. .I buf
  353. of length
  354. .I blen
  355. to an
  356. .IR mpint .
  357. If
  358. .I b
  359. is not
  360. .IR nil ,
  361. it refers to a preallocated
  362. .I mpint
  363. for the result.
  364. If
  365. .I b
  366. is
  367. .BR nil ,
  368. a new integer is allocated and returned as the result.
  369. .PP
  370. The integer conversions are:
  371. .TF Mptouv
  372. .TP
  373. .I mptoui
  374. .BR mpint -> "unsigned int"
  375. .TP
  376. .I uitomp
  377. .BR "unsigned int" -> mpint
  378. .TP
  379. .I mptoi
  380. .BR mpint -> "int"
  381. .TP
  382. .I itomp
  383. .BR "int" -> mpint
  384. .TP
  385. .I mptouv
  386. .BR mpint -> "unsigned vlong"
  387. .TP
  388. .I uvtomp
  389. .BR "unsigned vlong" -> mpint
  390. .TP
  391. .I mptov
  392. .BR mpint -> "vlong"
  393. .TP
  394. .I vtomp
  395. .BR "vlong" -> mpint
  396. .PD
  397. .PP
  398. When converting to the base integer types, if the integer is too large,
  399. the largest integer of the appropriate sign
  400. and size is returned.
  401. .PP
  402. The mathematical functions are:
  403. .TF mpmagadd
  404. .TP
  405. .I mpadd
  406. .BR "sum = b1 + b2" .
  407. .TP
  408. .I mpmagadd
  409. .BR "sum = abs(b1) + abs(b2)" .
  410. .TP
  411. .I mpsub
  412. .BR "diff = b1 - b2" .
  413. .TP
  414. .I mpmagsub
  415. .BR "diff = abs(b1) - abs(b2)" .
  416. .TP
  417. .I mpleft
  418. .BR "res = b<<shift" .
  419. .TP
  420. .I mpright
  421. .BR "res = b>>shift" .
  422. .TP
  423. .I mpmul
  424. .BR "prod = b1*b2" .
  425. .TP
  426. .I mpexp
  427. if
  428. .I m
  429. is nil,
  430. .BR "res = b**e" .
  431. Otherwise,
  432. .BR "res = b**e mod m" .
  433. .TP
  434. .I mpmod
  435. .BR "remainder = b % m" .
  436. .TP
  437. .I mpdiv
  438. .BR "quotient = dividend/divisor" .
  439. .BR "remainder = dividend % divisor" .
  440. .TP
  441. .I mpcmp
  442. returns -1, 0, or +1 as
  443. .I b1
  444. is less than, equal to, or greater than
  445. .IR b2 .
  446. .TP
  447. .I mpmagcmp
  448. the same as
  449. .I mpcmp
  450. but ignores the sign and just compares magnitudes.
  451. .PD
  452. .PP
  453. .I Mpextendedgcd
  454. computes the greatest common denominator,
  455. .IR d ,
  456. of
  457. .I a
  458. and
  459. .IR b .
  460. It also computes
  461. .I x
  462. and
  463. .I y
  464. such that
  465. .BR "a*x + b*y = d" .
  466. Both
  467. .I a
  468. and
  469. .I b
  470. are required to be positive.
  471. If called with negative arguments, it will
  472. return a gcd of 0.
  473. .PP
  474. .I Mpinverse
  475. computes the multiplicative inverse of
  476. .I b
  477. .B mod
  478. .IR m .
  479. .PP
  480. .I Mpsignif
  481. returns the bit offset of the left most 1 bit in
  482. .IR b .
  483. .I Mplowbits0
  484. returns the bit offset of the right most 1 bit.
  485. For example, for 0x14,
  486. .I mpsignif
  487. would return 4 and
  488. .I mplowbits0
  489. would return 2.
  490. .PP
  491. The remaining routines all work on arrays of
  492. .B mpdigit
  493. rather than
  494. .BR mpint 's.
  495. They are the basis of all the other routines. They are separated out
  496. to allow them to be rewritten in assembler for each architecture. There
  497. is also a portable C version for each one.
  498. .TF mpvecdigmuladd
  499. .TP
  500. .I mpdigdiv
  501. .BR "quotient = dividend[0:1] / divisor" .
  502. .TP
  503. .I mpvecadd
  504. .BR "sum[0:alen] = a[0:alen-1] + b[0:blen-1]" .
  505. We assume alen >= blen and that sum has room for alen+1 digits.
  506. .TP
  507. .I mpvecsub
  508. .BR "diff[0:alen-1] = a[0:alen-1] - b[0:blen-1]" .
  509. We assume that alen >= blen and that diff has room for alen digits.
  510. .TP
  511. .I mpvecdigmuladd
  512. .BR "p[0:n] += m * b[0:n-1]" .
  513. This multiplies a an array of digits times a scalar and adds it to another array.
  514. We assume p has room for n+1 digits.
  515. .TP
  516. .I mpvecdigmulsub
  517. .BR "p[0:n] -= m * b[0:n-1]" .
  518. This multiplies a an array of digits times a scalar and subtracts it fromo another array.
  519. We assume p has room for n+1 digits. It returns +1 is the result is positive and
  520. -1 if negative.
  521. .TP
  522. .I mpvecmul
  523. .BR "p[0:alen*blen] = a[0:alen-1] * b[0:blen-1]" .
  524. We assume that p has room for alen*blen+1 digits.
  525. .TP
  526. .I mpveccmp
  527. This returns -1, 0, or +1 as a - b is negative, 0, or positive.
  528. .PD
  529. .PP
  530. .IR mptwo ,
  531. .I mpone
  532. and
  533. .I mpzero
  534. are the constants 2, 1 and 0. These cannot be freed.
  535. .SS "Chinese remainder theorem
  536. .PP
  537. When computing in a non-prime modulus,
  538. .IR n,
  539. it is possible to perform the computations on the residues modulo the prime
  540. factors of
  541. .I n
  542. instead. Since these numbers are smaller, multiplication and exponentiation
  543. can be much faster.
  544. .PP
  545. .I Crtin
  546. computes the residues of
  547. .I x
  548. and returns them in a newly allocated structure:
  549. .EX
  550. typedef struct CRTres CRTres;
  551. {
  552. int n; // number of residues
  553. mpint *r[n]; // residues
  554. };
  555. .EE
  556. .PP
  557. .I Crtout
  558. takes a residue representation of a number and converts it back into
  559. the number. It also frees the residue structure.
  560. .PP
  561. .I Crepre
  562. saves a copy of the factors and precomputes the constants necessary
  563. for converting the residue form back into a number modulo
  564. the product of the factors. It returns a newly allocated structure
  565. containing values.
  566. .PP
  567. .I Crtprefree
  568. and
  569. .I crtresfree
  570. free
  571. .I CRTpre
  572. and
  573. .I CRTres
  574. structures respectively.
  575. .SH SOURCE
  576. .B /sys/src/libmp