1 /* mpz_gcdext(g, s, t, a, b) -- Set G to gcd(a, b), and S and T such that
4 Copyright 1991, 1993, 1994, 1995, 1996, 1997, 2000, 2001, 2005 Free Software
7 This file is part of the GNU MP Library.
9 The GNU MP Library is free software; you can redistribute it and/or modify
10 it under the terms of the GNU Lesser General Public License as published by
11 the Free Software Foundation; either version 3 of the License, or (at your
12 option) any later version.
14 The GNU MP Library is distributed in the hope that it will be useful, but
15 WITHOUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITY
16 or FITNESS FOR A PARTICULAR PURPOSE. See the GNU Lesser General Public
17 License for more details.
19 You should have received a copy of the GNU Lesser General Public License
20 along with the GNU MP Library. If not, see http://www.gnu.org/licenses/. */
22 #include <stdio.h> /* for NULL */
27 mpz_gcdext (mpz_ptr g
, mpz_ptr s
, mpz_ptr t
, mpz_srcptr a
, mpz_srcptr b
)
29 mp_size_t asize
, bsize
, usize
, vsize
;
32 mp_size_t gsize
, ssize
, tmp_ssize
;
33 mp_ptr gp
, sp
, tmp_gp
, tmp_sp
;
36 __mpz_struct stmp
, gtmp
;
41 /* mpn_gcdext requires that U >= V. Therefore, we often have to swap U and
42 V. This in turn leads to a lot of complications. The computed cofactor
43 will be the wrong one, so we have to fix that up at the end. */
45 asize
= ABS (SIZ (a
));
46 bsize
= ABS (SIZ (b
));
49 if (asize
> bsize
|| (asize
== bsize
&& mpn_cmp (ap
, bp
, asize
) > 0))
53 up
= TMP_ALLOC_LIMBS (usize
+ 1);
54 vp
= TMP_ALLOC_LIMBS (vsize
+ 1);
55 MPN_COPY (up
, ap
, usize
);
56 MPN_COPY (vp
, bp
, vsize
);
66 up
= TMP_ALLOC_LIMBS (usize
+ 1);
67 vp
= TMP_ALLOC_LIMBS (vsize
+ 1);
68 MPN_COPY (up
, bp
, usize
);
69 MPN_COPY (vp
, ap
, vsize
);
76 tmp_gp
= TMP_ALLOC_LIMBS (usize
+ 1);
77 tmp_sp
= TMP_ALLOC_LIMBS (usize
+ 1);
83 MPN_COPY (tmp_gp
, up
, usize
);
87 gsize
= mpn_gcdext (tmp_gp
, tmp_sp
, &tmp_ssize
, up
, usize
, vp
, vsize
);
88 ssize
= ABS (tmp_ssize
);
94 SIZ (&stmp
) = (tmp_ssize
^ SIZ (u
)) >= 0 ? ssize
: -ssize
;
103 MPZ_TMP_INIT (x
, ssize
+ usize
+ 1);
104 mpz_mul (x
, &stmp
, u
);
105 mpz_sub (x
, >mp
, x
);
106 mpz_tdiv_q (tt
, x
, v
);
112 if (ALLOC (ss
) < ssize
)
113 _mpz_realloc (ss
, ssize
);
115 MPN_COPY (sp
, tmp_sp
, ssize
);
116 SIZ (ss
) = SIZ (&stmp
);
119 if (ALLOC (g
) < gsize
)
120 _mpz_realloc (g
, gsize
);
122 MPN_COPY (gp
, tmp_gp
, gsize
);