#include #include "gmp_additional_functions.h" void mpz_inits (mpz_ptr x, ...) { va_list ap; va_start (ap, x); do { mpz_init(x); x = va_arg (ap, mpz_ptr); }while (x != NULL); va_end (ap); } void mpz_clears (mpz_ptr x, ...) { va_list ap; va_start (ap, x); do { mpz_clear(x); x = va_arg (ap, mpz_ptr); }while (x != NULL); va_end (ap); } void mpq_inits (mpq_ptr x, ...) { va_list ap; va_start (ap, x); do { mpq_init(x); x = va_arg (ap, mpq_ptr); }while (x != NULL); va_end (ap); } void mpq_clears (mpq_ptr x, ...) { va_list ap; va_start (ap, x); do { mpq_clear(x); x = va_arg (ap, mpq_ptr); }while (x != NULL); va_end (ap); } // https://github.com/python/cpython/blob/3.10/Lib/fractions.py#L249 int limit_denominator (mpq_t fraction, mpz_t max_denominator) { mpz_t p0, q0, p1, q1, q2; mpz_t n, d; mpz_t a, temp; __mpz_struct *numerator = mpq_numref (fraction); __mpz_struct *denominator = mpq_denref (fraction); if (mpz_sgn(max_denominator) <= 0) return(1); if (mpz_cmp(denominator,max_denominator) <= 0) return(0); mpz_inits (p0, q0, p1, q1, q2, n, d, a, temp,NULL); mpz_set_ui(p0,0); mpz_set_ui(q0,1); mpz_set_ui(p1,1); mpz_set_ui(q1,0); // n, d = self._numerator, self._denominator mpz_set(n,numerator); mpz_set(d,denominator); while(1) { // a = n//d // if (!mpz_sgn(d)) break; //necesery only if input is not canonicalized mpz_fdiv_q (a,n,d); // q2 = q0+a*q1 mpz_set (q2, q0); mpz_addmul (q0, a, q1); // if q2 > max_denominator: // break if (mpz_cmp(q0,max_denominator) > 0) break; // p0, q0, p1, q1 = p1, q1, p0+a*p1, q2 mpz_addmul (p0, a, p1); mpz_swap (p0, p1); mpz_swap (q0,q1); // n, d = d, n-a*d mpz_submul (n,a,d); mpz_swap (n, d); } mpz_t k, b1d, b1n; mpz_inits (k, b1d, b1n, NULL); // k = (max_denominator-q0)//q1 mpz_sub(temp,max_denominator,q2); mpz_tdiv_q (k,temp,q1); // bound1 = Fraction(p0+k*p1, q0+k*q1) // bound2 = Fraction(p1, q1) mpz_set(b1d,q2); mpz_addmul(b1d,k,q1); mpz_set(b1n,p0); mpz_addmul(b1n,k,p1); mpz_t b1_cmp, b2_cmp, self_cmp, b1_d, b2_d, r_n, r_d; mpz_inits (b1_cmp, b2_cmp, self_cmp, b1_d, b2_d, r_n, r_d, NULL); // if abs(bound2 - self) <= abs(bound1-self): // return bound2 // else: // return bound1 mpz_mul(temp,q1,denominator); mpz_mul(b1_cmp,temp,b1n); mpz_mul(temp,b1d,denominator); mpz_mul(b2_cmp,temp,p1); mpz_mul(temp,q1,b1d); mpz_mul(self_cmp,temp,numerator); mpz_sub(b1_d,b1_cmp,self_cmp); mpz_sub(b2_d,b2_cmp,self_cmp); if (mpz_cmpabs(b2_d,b1_d) <= 0) { mpz_set(numerator,p1); mpz_set(denominator,q1); } else { mpz_set(numerator,b1n); mpz_set(denominator,b1d); } mpz_clears (p0, q0, p1, q1, q2, n, d, a, temp,NULL); mpz_clears (k, b1d, b1n, NULL); mpz_clears (b1_cmp, b2_cmp, self_cmp, b1_d, b2_d, r_n, r_d, NULL); return(0); } // void fprint_mpq (mpq_t value, unsigned int size, FILE *fd) // { // mpz_t temp, temp2; // mpz_inits (temp, temp2, NULL); // __mpz_struct *value_n = mpq_numref (value); // __mpz_struct *value_d = mpq_denref (value); // if (mpq_sgn (value) < 0) // { // fprintf (fd, "-"); // } // mpz_tdiv_q (temp, value_n, value_d); // mpz_abs (temp, temp); // mpz_out_str (fd, 10, temp); // mpz_tdiv_r (temp2, value_n, value_d); // mpz_abs (temp2, temp2); // mpz_ui_pow_ui (temp, 10, size); // mpz_mul (temp, temp2, temp); // mpz_tdiv_q_ui (temp2, value_d, 2); // mpz_add (temp, temp, temp2); // mpz_tdiv_q (temp, temp, value_d); // gmp_fprintf (fd, ".%0*Zd", size, temp); // mpz_clears (temp, temp2, NULL); // } // void print_mpq (mpq_t value, unsigned int size) // { // fprint_mpq (value, size, stdout); // } int mpq_from_str_internal (mpq_t value, char* str) { __mpz_struct *value_n = mpq_numref (value); __mpz_struct *value_d = mpq_denref (value); unsigned int i; char *buff; unsigned int size; unsigned int start; char separator; for (i = 0; (str[i] >= '0' && str[i] <= '9'); i++); separator = str[i]; size = i; if (separator == '\0') //integer { if (size == 0) return (1); mpz_set_str (value_n, str, 10); mpz_set_ui (value_d, 1); return (0); } if (!(separator == '.' || separator == ',' || separator == '/')) return (1); if (size == 0) return (1); buff = (char*) malloc(size + 1); memcpy (buff, str, size); buff[size] = '\0'; mpz_set_str (value_n, buff, 10); free (buff); start = size + 1; for (i = start; (str[i] >= '0' && str[i] <= '9'); i++); if (str[i] != '\0') return (1); size = i - start; if (size == 0) return (1); if (separator == '/') //fraction { mpz_set_str (value_d, str + start, 10); return (0); } mpz_t temp; mpz_init (temp); mpz_set_str (temp, str + start, 10); mpz_ui_pow_ui (value_d, 10, size); mpz_mul (value_n, value_n, value_d); mpz_add (value_n, value_n, temp); mpz_clear (temp); return (0); } int mpq_from_str_internal2 (mpq_t value, char* str) { int res = mpq_from_str_internal (value, str); if (res) mpq_set_ui (value, 0, 1); mpq_canonicalize (value); return (res); } int mpq_set_str_plus (mpq_t value, char* str) { int size = strlen (str); char n = 0; if (size == 0) return (1); switch (str[size - 1]) { case'p': n = -12; break; case'n': n = -9; break; case'u': n = -6; break; case'm': n = -3; break; case'k': case'K': n = 3; break; case'M': n = 6; break; case'G': n = 9; break; case'T': n = 12; break; } if (n) { // printf("aaa\n"); char *str2 = (char*)malloc(size); strncpy(str2, str, size - 1); str2[size - 1] = '\0'; // Null-terminate the string int ret = mpq_from_str_internal2 (value, str2); free(str2); mpq_t temp; mpq_init (temp); mpq_set_ui (temp, 1, 1); if (n > 0) { mpz_ui_pow_ui (mpq_numref (temp), 10, n); } else { mpz_ui_pow_ui (mpq_denref (temp), 10, -n); } mpq_mul (value, value, temp); mpq_clear (temp); return (ret); } else { // printf("bbb%d\n",str[size - 1]); return (mpq_from_str_internal2 (value, str)); } } // int limit_denominator_alt (mpq_t fraction, mpz_t max_denominator, mpq_t alt_fraction) // { // mpz_t p0, q0, p1, q1, q2; // mpz_t n, d; // mpz_t a, temp; // // __mpz_struct *numerator = mpq_numref (fraction); // __mpz_struct *denominator = mpq_denref (fraction); // __mpz_struct *alt_numerator = mpq_numref (alt_fraction); // __mpz_struct *alt_denominator = mpq_denref (alt_fraction); // // if (mpz_sgn(max_denominator) <= 0) return(1); // if (mpz_cmp(denominator,max_denominator) <= 0) // { // mpq_set (alt_fraction, fraction); // return(0); // } // // mpz_inits (p0, q0, p1, q1, q2, n, d, a, temp,NULL); // // mpz_set_ui(p0,0); // mpz_set_ui(q0,1); // mpz_set_ui(p1,1); // mpz_set_ui(q1,0); // // // n, d = self._numerator, self._denominator // // mpz_set(n,numerator); // mpz_set(d,denominator); // // while(1) // { // // a = n//d // // if (!mpz_sgn(d)) break; //necesery only if input is not canonicalized // mpz_fdiv_q (a,n,d); // // // q2 = q0+a*q1 // mpz_set (q2, q0); // mpz_addmul (q0, a, q1); // // // if q2 > max_denominator: // // break // if (mpz_cmp(q0,max_denominator) > 0) break; // // // p0, q0, p1, q1 = p1, q1, p0+a*p1, q2 // mpz_addmul (p0, a, p1); // mpz_swap (p0, p1); // // mpz_swap (q0,q1); // // // n, d = d, n-a*d // mpz_submul (n,a,d); // mpz_swap (n, d); // } // // mpz_t k, b1d, b1n; // // mpz_inits (k, b1d, b1n, NULL); // // // k = (max_denominator-q0)//q1 // mpz_sub(temp,max_denominator,q2); // mpz_tdiv_q (k,temp,q1); // // // bound1 = Fraction(p0+k*p1, q0+k*q1) // // bound2 = Fraction(p1, q1) // mpz_set(b1d,q2); // mpz_addmul(b1d,k,q1); // // mpz_set(b1n,p0); // mpz_addmul(b1n,k,p1); // // mpz_t b1_cmp, b2_cmp, self_cmp, b1_d, b2_d, r_n, r_d; // mpz_inits (b1_cmp, b2_cmp, self_cmp, b1_d, b2_d, r_n, r_d, NULL); // // // if abs(bound2 - self) <= abs(bound1-self): // // return bound2 // // else: // // return bound1 // // mpz_mul(temp,q1,denominator); // mpz_mul(b1_cmp,temp,b1n); // // mpz_mul(temp,b1d,denominator); // mpz_mul(b2_cmp,temp,p1); // // mpz_mul(temp,q1,b1d); // mpz_mul(self_cmp,temp,numerator); // // mpz_sub(b1_d,b1_cmp,self_cmp); // mpz_sub(b2_d,b2_cmp,self_cmp); // // if (mpz_cmpabs(b2_d,b1_d) <= 0) // { // mpz_set(numerator,p1); // mpz_set(denominator,q1); // // mpz_set(alt_numerator,b1n); // mpz_set(alt_denominator,b1d); // } // else // { // mpz_set(numerator,b1n); // mpz_set(denominator,b1d); // // mpz_set(alt_numerator,p1); // mpz_set(alt_denominator,q1); // } // // mpz_clears (p0, q0, p1, q1, q2, n, d, a, temp,NULL); // mpz_clears (k, b1d, b1n, NULL); // mpz_clears (b1_cmp, b2_cmp, self_cmp, b1_d, b2_d, r_n, r_d, NULL); // // return(0); // }