1212#include "flint.h"
1313#include "gr.h"
1414#include "gr_ore_poly.h"
15+ #include "gr_vec.h"
1516
1617// Returns the unique pair (Q, R) such that U = QV + R and ord(R) < ord(V)
1718int _gr_ore_poly_divrem (gr_ptr Q , gr_ptr R , gr_srcptr U , slong lenU , gr_srcptr V , slong lenV , gr_ore_poly_ctx_t ctx )
@@ -21,78 +22,61 @@ int _gr_ore_poly_divrem(gr_ptr Q, gr_ptr R, gr_srcptr U, slong lenU, gr_srcptr V
2122 slong lenQ , lenR , ordV , ordR ;
2223 int status = GR_SUCCESS ;
2324
24- if (GR_ORE_POLY_CTX (ctx )-> sigma_delta == NULL )
25- return GR_UNABLE ;
26-
2725 lenQ = lenU - lenV + 1 ;
2826 lenR = lenU ;
2927 ordV = lenV - 1 ;
3028 ordR = lenR - 1 ;
3129
3230 // Set Q to 0
33- for (slong k = 0 ; k < lenQ ; k ++ )
34- status |= gr_zero (GR_ENTRY (Q , k , el_size ), cctx );
31+ status |= _gr_vec_zero (Q , lenQ , cctx );
3532
3633 // Set R to U
37- for (slong k = 0 ; k < lenU ; k ++ )
38- status |= gr_set (GR_ENTRY (R , k , el_size ), GR_ENTRY (U , k , el_size ), cctx );
34+ status |= _gr_vec_set (R , U , lenU , cctx );
3935
40- gr_ptr lcR , lcV , denominator , c ;
41- GR_TMP_INIT (lcR , cctx );
42- GR_TMP_INIT (lcV , cctx );
43- GR_TMP_INIT (denominator , cctx );
36+ gr_ptr c ;
4437 GR_TMP_INIT (c , cctx );
4538
46- gr_ptr A = flint_malloc (lenQ * el_size );
47- gr_ptr B = flint_malloc (lenU * el_size );
48- _gr_vec_init (A , lenQ , cctx );
49- _gr_vec_init (B , lenU , cctx );
39+ gr_ptr A , B ;
40+ GR_TMP_INIT_VEC (A , lenQ , cctx );
41+ GR_TMP_INIT_VEC (B , lenU , cctx );
42+
43+ gr_ptr sigma_pows ;
44+ GR_TMP_INIT_VEC (sigma_pows , lenQ , cctx );
45+ status |= gr_set (GR_ENTRY (sigma_pows , 0 , el_size ), GR_ENTRY (V , ordV , el_size ), cctx );
46+ for (slong i = 1 ; i < lenQ ; i ++ )
47+ status |= gr_ore_poly_sigma (GR_ENTRY (sigma_pows , i , el_size ), GR_ENTRY (sigma_pows , i - 1 , el_size ), ctx );
5048
5149 while (ordR > ordV )
5250 {
5351 slong k = ordR - ordV ;
5452
55- status |= gr_set (lcR , GR_ENTRY (R , ordR , el_size ), cctx );
56- status |= gr_set (lcV , GR_ENTRY (V , ordV , el_size ), cctx );
57-
58- // Compute denominator = sigma ^ k (lc(V))
59- status |= gr_set (denominator , lcV , cctx );
60- for (slong i = 0 ; i < k ; i ++ )
61- status |= gr_ore_poly_sigma (denominator , denominator , ctx );
62-
63- // c = lc(R) / denominator
64- status |= gr_div (c , lcR , denominator , cctx );
53+ // c = lc(R) / sigma ^ k (lc(V)) using sigma_pows computed above
54+ status |= gr_div (c , GR_ENTRY (R , ordR , el_size ), GR_ENTRY (sigma_pows , k , el_size ), cctx );
6555
6656 // R -= c * x^k * V. We compute A = c * x^k, then B = A * V
6757 // A = c * x^k, so A[k] = c, rest 0
6858 slong lenA = k + 1 ;
69- for (slong i = 0 ; i < lenA ; i ++ )
70- status |= gr_zero (GR_ENTRY (A , i , el_size ), cctx );
59+ status |= _gr_vec_zero (A , lenA , cctx );
7160 status |= gr_set (GR_ENTRY (A , k , el_size ), c , cctx );
7261
7362 // B = A * V
7463 status |= _gr_ore_poly_mul (B , A , lenA , V , lenV , ctx );
7564
7665 // R -= B
7766 slong lenB = lenA + lenV - 1 ;
78- for (slong i = 0 ; i < lenB ; i ++ )
79- status |= gr_sub (GR_ENTRY (R , i , el_size ), GR_ENTRY (R , i , el_size ), GR_ENTRY (B , i , el_size ), cctx );
80-
67+ status |= _gr_vec_sub (R , R , B , lenB , cctx );
68+
8169 // Q += c * x^k, so Q[k] += c
8270 status |= gr_add (GR_ENTRY (Q , k , el_size ), GR_ENTRY (Q , k , el_size ), c , cctx );
8371
8472 ordR -- ;
8573 }
8674
87- gr_clear (lcR , cctx );
88- gr_clear (lcV , cctx );
89- gr_clear (denominator , cctx );
90- gr_clear (c , cctx );
75+ GR_TMP_CLEAR_VEC (sigma_pows , lenQ , cctx );
76+ GR_TMP_CLEAR (c , cctx );
9177
92- _gr_vec_clear (A , lenQ , cctx );
93- _gr_vec_clear (B , lenU , cctx );
94- flint_free (A );
95- flint_free (B );
78+ GR_TMP_CLEAR_VEC (A , lenQ , cctx );
79+ GR_TMP_CLEAR_VEC (B , lenU , cctx );
9680
9781 return status ;
9882}
@@ -115,19 +99,14 @@ int gr_ore_poly_divrem(gr_ore_poly_t Q, gr_ore_poly_t R, const gr_ore_poly_t U,
11599
116100 if (lenU < lenV )
117101 {
118- gr_ore_poly_t tU ;
119- // take care of aliasing case with temp
120- gr_ore_poly_init (tU , ctx );
121- status |= gr_ore_poly_set (tU , U , ctx );
102+ status |= gr_ore_poly_set (R , U , ctx );
122103 status |= gr_ore_poly_zero (Q , ctx );
123- status |= gr_ore_poly_set (R , tU , ctx );
124- gr_ore_poly_clear (tU , ctx );
125104 return status ;
126105 }
127106
128107 slong lenQ = lenU - lenV + 1 ;
129108
130- if (Q == U || Q == V || R == U || R == V ) // treat aliasing case separately
109+ if (Q == U || Q == V || R == U || R == V || Q == R ) // treat aliasing case separately
131110 {
132111 gr_ore_poly_t tQ , tR ;
133112 gr_ore_poly_init (tQ , ctx );
@@ -158,3 +137,23 @@ int gr_ore_poly_divrem(gr_ore_poly_t Q, gr_ore_poly_t R, const gr_ore_poly_t U,
158137 }
159138 return status ;
160139}
140+
141+ int gr_ore_poly_div (gr_ore_poly_t Q , const gr_ore_poly_t U , gr_ore_poly_t V , gr_ore_poly_ctx_t ctx )
142+ {
143+ gr_ore_poly_t R ;
144+ int status ;
145+ gr_ore_poly_init (R , ctx );
146+ status = gr_ore_poly_divrem (Q , R , U , V , ctx );
147+ gr_ore_poly_clear (R , ctx );
148+ return status ;
149+ }
150+
151+ int gr_ore_poly_rem (gr_ore_poly_t R , const gr_ore_poly_t U , gr_ore_poly_t V , gr_ore_poly_ctx_t ctx )
152+ {
153+ gr_ore_poly_t Q ;
154+ int status ;
155+ gr_ore_poly_init (Q , ctx );
156+ status = gr_ore_poly_divrem (Q , R , U , V , ctx );
157+ gr_ore_poly_clear (Q , ctx );
158+ return status ;
159+ }
0 commit comments