diff options
| author | Kaleb Keithley <kaleb@freedesktop.org> | 2003-11-14 16:49:22 +0000 |
|---|---|---|
| committer | Kaleb Keithley <kaleb@freedesktop.org> | 2003-11-14 16:49:22 +0000 |
| commit | 0a193e032ba1ecf3f003e027e833dc9d274cb740 (patch) | |
| tree | a1dcc00cb7f5d26e437e05e658c38fc323fe919d /lisp/mathimp.c | |
Initial revision
Diffstat (limited to 'lisp/mathimp.c')
| -rw-r--r-- | lisp/mathimp.c | 5225 |
1 files changed, 5225 insertions, 0 deletions
diff --git a/lisp/mathimp.c b/lisp/mathimp.c new file mode 100644 index 0000000..ccda576 --- /dev/null +++ b/lisp/mathimp.c @@ -0,0 +1,5225 @@ +/* + * Copyright (c) 2002 by The XFree86 Project, Inc. + * + * Permission is hereby granted, free of charge, to any person obtaining a + * copy of this software and associated documentation files (the "Software"), + * to deal in the Software without restriction, including without limitation + * the rights to use, copy, modify, merge, publish, distribute, sublicense, + * and/or sell copies of the Software, and to permit persons to whom the + * Software is furnished to do so, subject to the following conditions: + * + * The above copyright notice and this permission notice shall be included in + * all copies or substantial portions of the Software. + * + * THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS OR + * IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY, + * FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL + * THE XFREE86 PROJECT BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY, + * WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING FROM, OUT OF + * OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE + * SOFTWARE. + * + * Except as contained in this notice, the name of the XFree86 Project shall + * not be used in advertising or otherwise to promote the sale, use or other + * dealings in this Software without prior written authorization from the + * XFree86 Project. + * + * Author: Paulo César Pereira de Andrade + */ + +/* $XFree86: xc/programs/xedit/lisp/mathimp.c,v 1.14 2003/01/30 02:46:25 paulo Exp $ */ + + +/* + * Defines + */ +#ifdef __GNUC__ +#define CONST __attribute__ ((__const__)) +#else +#define CONST /**/ +#endif + +/* mask for checking overflow on long operations */ +#ifdef LONG64 +#define FI_MASK 0x4000000000000000L +#define LONGSBITS 63 +#else +#define FI_MASK 0x40000000L +#define LONGSBITS 31 +#endif + +#define N_FIXNUM 1 +#define N_BIGNUM 2 +#define N_FLONUM 3 +#define N_FIXRATIO 4 +#define N_BIGRATIO 5 + +#define NOP_ADD 1 +#define NOP_SUB 2 +#define NOP_MUL 3 +#define NOP_DIV 4 + +#define NDIVIDE_CEIL 1 +#define NDIVIDE_FLOOR 2 +#define NDIVIDE_ROUND 3 +#define NDIVIDE_TRUNC 4 + +/* real part from number */ +#define NREAL(num) &((num)->real) +#define NRTYPE(num) (num)->real.type +#define NRFI(num) (num)->real.data.fixnum +#define NRBI(num) (num)->real.data.bignum +#define NRFF(num) (num)->real.data.flonum +#define NRFRN(Num) (Num)->real.data.fixratio.num +#define NRFRD(num) (num)->real.data.fixratio.den +#define NRBR(num) (num)->real.data.bigratio +#define NRBRN(num) mpr_num(NRBR(num)) +#define NRBRD(num) mpr_den(NRBR(num)) + +#define NRCLEAR_BI(num) mpi_clear(NRBI(num)); XFREE(NRBI(num)) +#define NRCLEAR_BR(num) mpr_clear(NRBR(num)); XFREE(NRBR(num)) + +/* imag part from number */ +#define NIMAG(num) &((num)->imag) +#define NITYPE(num) (num)->imag.type +#define NIFI(num) (num)->imag.data.fixnum +#define NIBI(num) (num)->imag.data.bignum +#define NIFF(num) (num)->imag.data.flonum +#define NIFRN(Num) (Num)->imag.data.fixratio.num +#define NIFRD(num) (num)->imag.data.fixratio.den +#define NIBR(num) (num)->imag.data.bigratio +#define NIBRN(obj) mpr_num(NIBR(obj)) +#define NIBRD(obj) mpr_den(NIBR(obj)) + +/* real number fields */ +#define RTYPE(real) (real)->type +#define RFI(real) (real)->data.fixnum +#define RBI(real) (real)->data.bignum +#define RFF(real) (real)->data.flonum +#define RFRN(real) (real)->data.fixratio.num +#define RFRD(real) (real)->data.fixratio.den +#define RBR(real) (real)->data.bigratio +#define RBRN(real) mpr_num(RBR(real)) +#define RBRD(real) mpr_den(RBR(real)) + +#define RINTEGERP(real) \ + (RTYPE(real) == N_FIXNUM || RTYPE(real) == N_BIGNUM) + +#define RCLEAR_BI(real) mpi_clear(RBI(real)); XFREE(RBI(real)) +#define RCLEAR_BR(real) mpr_clear(RBR(real)); XFREE(RBR(real)) + +/* numeric value from lisp object */ +#define OFI(object) FIXNUM_VALUE(object) +#define OII(object) INT_VALUE(object) +#define OBI(object) (object)->data.mp.integer +#define ODF(object) DFLOAT_VALUE(object) +#define OFRN(object) (object)->data.ratio.numerator +#define OFRD(object) (object)->data.ratio.denominator +#define OBR(object) (object)->data.mp.ratio +#define OBRN(object) mpr_num(OBR(object)) +#define OBRD(object) mpr_den(OBR(object)) +#define OCXR(object) (object)->data.complex.real +#define OCXI(object) (object)->data.complex.imag + +#define XALLOC(type) LispMalloc(sizeof(type)) +#define XFREE(ptr) LispFree(ptr) + + +/* + * Types + */ +typedef struct _n_real { + char type; + union { + long fixnum; + mpi *bignum; + double flonum; + struct { + long num; + long den; + } fixratio; + mpr *bigratio; + } data; +} n_real; + +typedef struct _n_number { + char complex; + n_real real; + n_real imag; +} n_number; + + +/* + * Prototypes + */ +static void number_init(void); +static LispObj *number_pi(void); + +static void set_real_real(n_real*, n_real*); +static void set_real_object(n_real*, LispObj*); +static void set_number_object(n_number*, LispObj*); +static void clear_real(n_real*); +static void clear_number(n_number*); + +static LispObj *make_real_object(n_real*); +static LispObj *make_number_object(n_number*); + +static void fatal_error(int); +static void fatal_object_error(LispObj*, int); +static void fatal_builtin_object_error(LispBuiltin*, LispObj*, int); + +static double bi_getd(mpi*); +static double br_getd(mpr*); + +/* add */ +static void add_real_object(n_real*, LispObj*); +static void add_number_object(n_number*, LispObj*); + +/* sub */ +static void sub_real_object(n_real*, LispObj*); +static void sub_number_object(n_number*, LispObj*); + +/* mul */ +static void mul_real_object(n_real*, LispObj*); +static void mul_number_object(n_number*, LispObj*); + +/* div */ +static void div_real_object(n_real*, LispObj*); +static void div_number_object(n_number*, LispObj*); + +/* compare */ +static int cmp_real_real(n_real*, n_real*); +static int cmp_real_object(n_real*, LispObj*); +#if 0 /* not used */ +static int cmp_number_object(n_number*, LispObj*); +#endif +static int cmp_object_object(LispObj*, LispObj*, int); + +/* fixnum */ +static INLINE int fi_fi_add_overflow(long, long) CONST; +static INLINE int fi_fi_sub_overflow(long, long) CONST; +static INLINE int fi_fi_mul_overflow(long, long) CONST; + +/* bignum */ +static void rbi_canonicalize(n_real*); + +/* ratio */ +static void rfr_canonicalize(n_real*); +static void rbr_canonicalize(n_real*); + +/* complex */ +static void ncx_canonicalize(n_number*); + +/* abs */ +static void abs_real(n_real*); +static void abs_number(n_number*); +static void nabs_cx(n_number*); +static INLINE void rabs_fi(n_real*); +static INLINE void rabs_bi(n_real*); +static INLINE void rabs_ff(n_real*); +static INLINE void rabs_fr(n_real*); +static INLINE void rabs_br(n_real*); + +/* neg */ +static void neg_real(n_real*); +static void neg_number(n_number*); +static void rneg_fi(n_real*); +static INLINE void rneg_bi(n_real*); +static INLINE void rneg_ff(n_real*); +static INLINE void rneg_fr(n_real*); +static INLINE void rneg_br(n_real*); + +/* sqrt */ +static void sqrt_real(n_real*); +static void sqrt_number(n_number*); +static void rsqrt_xi(n_real*); +static void rsqrt_xr(n_real*); +static void rsqrt_ff(n_real*); +static void nsqrt_cx(n_number*); +static void nsqrt_xi(n_number*); +static void nsqrt_ff(n_number*); +static void nsqrt_xr(n_number*); + +/* mod */ +static void mod_real_real(n_real*, n_real*); +static void mod_real_object(n_real*, LispObj*); +static void rmod_fi_fi(n_real*, long); +static void rmod_fi_bi(n_real*, mpi*); +static void rmod_bi_fi(n_real*, long); +static void rmod_bi_bi(n_real*, mpi*); + +/* rem */ +static void rem_real_object(n_real*, LispObj*); +static void rrem_fi_fi(n_real*, long); +static void rrem_fi_bi(n_real*, mpi*); +static void rrem_bi_fi(n_real*, long); +static void rrem_bi_bi(n_real*, mpi*); + +/* gcd */ +static void gcd_real_object(n_real*, LispObj*); + +/* and */ +static void and_real_object(n_real*, LispObj*); + +/* eqv */ +static void eqv_real_object(n_real*, LispObj*); + +/* ior */ +static void ior_real_object(n_real*, LispObj*); + +/* not */ +static void not_real(n_real*); + +/* xor */ +static void xor_real_object(n_real*, LispObj*); + +/* divide */ +static void divide_number_object(n_number*, LispObj*, int, int); +static void ndivide_xi_xi(n_number*, LispObj*, int, int); +static void ndivide_flonum(n_number*, double, double, int, int); +static void ndivide_xi_xr(n_number*, LispObj*, int, int); +static void ndivide_xr_xi(n_number*, LispObj*, int, int); +static void ndivide_xr_xr(n_number*, LispObj*, int, int); + +/* real complex */ +static void nadd_re_cx(n_number*, LispObj*); +static void nsub_re_cx(n_number*, LispObj*); +static void nmul_re_cx(n_number*, LispObj*); +static void ndiv_re_cx(n_number*, LispObj*); + +/* complex real */ +static void nadd_cx_re(n_number*, LispObj*); +static void nsub_cx_re(n_number*, LispObj*); +static void nmul_cx_re(n_number*, LispObj*); +static void ndiv_cx_re(n_number*, LispObj*); + +/* complex complex */ +static void nadd_cx_cx(n_number*, LispObj*); +static void nsub_cx_cx(n_number*, LispObj*); +static void nmul_cx_cx(n_number*, LispObj*); +static void ndiv_cx_cx(n_number*, LispObj*); +static int cmp_cx_cx(LispObj*, LispObj*); + +/* flonum flonum */ +static void radd_flonum(n_real*, double, double); +static void rsub_flonum(n_real*, double, double); +static void rmul_flonum(n_real*, double, double); +static void rdiv_flonum(n_real*, double, double); +static int cmp_flonum(double, double); + +/* fixnum fixnum */ +static void rop_fi_fi_bi(n_real*, long, int); +static INLINE void radd_fi_fi(n_real*, long); +static INLINE void rsub_fi_fi(n_real*, long); +static INLINE void rmul_fi_fi(n_real*, long); +static INLINE void rdiv_fi_fi(n_real*, long); +static INLINE int cmp_fi_fi(long, long); +static void ndivide_fi_fi(n_number*, long, int, int); + +/* fixnum bignum */ +static void rop_fi_bi_xi(n_real*, mpi*, int); +static INLINE void radd_fi_bi(n_real*, mpi*); +static INLINE void rsub_fi_bi(n_real*, mpi*); +static INLINE void rmul_fi_bi(n_real*, mpi*); +static void rdiv_fi_bi(n_real*, mpi*); +static INLINE int cmp_fi_bi(long, mpi*); + +/* fixnum fixratio */ +static void rop_fi_fr_as_xr(n_real*, long, long, int); +static void rop_fi_fr_md_xr(n_real*, long, long, int); +static INLINE void radd_fi_fr(n_real*, long, long); +static INLINE void rsub_fi_fr(n_real*, long, long); +static INLINE void rmul_fi_fr(n_real*, long, long); +static INLINE void rdiv_fi_fr(n_real*, long, long); +static INLINE int cmp_fi_fr(long, long, long); + +/* fixnum bigratio */ +static void rop_fi_br_as_xr(n_real*, mpr*, int); +static void rop_fi_br_md_xr(n_real*, mpr*, int); +static INLINE void radd_fi_br(n_real*, mpr*); +static INLINE void rsub_fi_br(n_real*, mpr*); +static INLINE void rmul_fi_br(n_real*, mpr*); +static INLINE void rdiv_fi_br(n_real*, mpr*); +static INLINE int cmp_fi_br(long, mpr*); + +/* bignum fixnum */ +static INLINE void radd_bi_fi(n_real*, long); +static INLINE void rsub_bi_fi(n_real*, long); +static INLINE void rmul_bi_fi(n_real*, long); +static void rdiv_bi_fi(n_real*, long); +static INLINE int cmp_bi_fi(mpi*, long); + +/* bignum bignum */ +static INLINE void radd_bi_bi(n_real*, mpi*); +static INLINE void rsub_bi_bi(n_real*, mpi*); +static INLINE void rmul_bi_bi(n_real*, mpi*); +static void rdiv_bi_bi(n_real*, mpi*); +static INLINE int cmp_bi_bi(mpi*, mpi*); + +/* bignum fixratio */ +static void rop_bi_fr_as_xr(n_real*, long, long, int); +static void rop_bi_fr_md_xr(n_real*, long, long, int); +static INLINE void radd_bi_fr(n_real*, long, long); +static INLINE void rsub_bi_fr(n_real*, long, long); +static INLINE void rmul_bi_fr(n_real*, long, long); +static INLINE void rdiv_bi_fr(n_real*, long, long); +static int cmp_bi_fr(mpi*, long, long); + +/* bignum bigratio */ +static void rop_bi_br_as_xr(n_real*, mpr*, int); +static void rop_bi_br_md_xr(n_real*, mpr*, int); +static INLINE void radd_bi_br(n_real*, mpr*); +static INLINE void rsub_bi_br(n_real*, mpr*); +static INLINE void rmul_bi_br(n_real*, mpr*); +static INLINE void rdiv_bi_br(n_real*, mpr*); +static int cmp_bi_br(mpi*, mpr*); + +/* fixratio fixnum */ +static void rop_fr_fi_as_xr(n_real*, long, int); +static void rop_fr_fi_md_xr(n_real*, long, int); +static INLINE void radd_fr_fi(n_real*, long); +static INLINE void rsub_fr_fi(n_real*, long); +static INLINE void rmul_fr_fi(n_real*, long); +static INLINE void rdiv_fr_fi(n_real*, long); +static INLINE int cmp_fr_fi(long, long, long); + +/* fixratio bignum */ +static void rop_fr_bi_as_xr(n_real*, mpi*, int); +static void rop_fr_bi_md_xr(n_real*, mpi*, int); +static INLINE void radd_fr_bi(n_real*, mpi*); +static INLINE void rsub_fr_bi(n_real*, mpi*); +static INLINE void rmul_fr_bi(n_real*, mpi*); +static INLINE void rdiv_fr_bi(n_real*, mpi*); +static int cmp_fr_bi(long, long, mpi*); + +/* fixratio fixratio */ +static void rop_fr_fr_as_xr(n_real*, long, long, int); +static void rop_fr_fr_md_xr(n_real*, long, long, int); +static INLINE void radd_fr_fr(n_real*, long, long); +static INLINE void rsub_fr_fr(n_real*, long, long); +static INLINE void rmul_fr_fr(n_real*, long, long); +static INLINE void rdiv_fr_fr(n_real*, long, long); +static INLINE int cmp_fr_fr(long, long, long, long); + +/* fixratio bigratio */ +static void rop_fr_br_asmd_xr(n_real*, mpr*, int); +static INLINE void radd_fr_br(n_real*, mpr*); +static INLINE void rsub_fr_br(n_real*, mpr*); +static INLINE void rmul_fr_br(n_real*, mpr*); +static INLINE void rdiv_fr_br(n_real*, mpr*); +static int cmp_fr_br(long, long, mpr*); + +/* bigratio fixnum */ +static void rop_br_fi_asmd_xr(n_real*, long, int); +static INLINE void radd_br_fi(n_real*, long); +static INLINE void rsub_br_fi(n_real*, long); +static INLINE void rmul_br_fi(n_real*, long); +static INLINE void rdiv_br_fi(n_real*, long); +static int cmp_br_fi(mpr*, long); + +/* bigratio bignum */ +static void rop_br_bi_as_xr(n_real*, mpi*, int); +static INLINE void radd_br_bi(n_real*, mpi*); +static INLINE void rsub_br_bi(n_real*, mpi*); +static INLINE void rmul_br_bi(n_real*, mpi*); +static INLINE void rdiv_br_bi(n_real*, mpi*); +static int cmp_br_bi(mpr*, mpi*); + +/* bigratio fixratio */ +static void rop_br_fr_asmd_xr(n_real*, long, long, int); +static INLINE void radd_br_fr(n_real*, long, long); +static INLINE void rsub_br_fr(n_real*, long, long); +static INLINE void rmul_br_fr(n_real*, long, long); +static INLINE void rdiv_br_fr(n_real*, long, long); +static int cmp_br_fr(mpr*, long, long); + +/* bigratio bigratio */ +static INLINE void radd_br_br(n_real*, mpr*); +static INLINE void rsub_br_br(n_real*, mpr*); +static INLINE void rmul_br_br(n_real*, mpr*); +static INLINE void rdiv_br_br(n_real*, mpr*); +static INLINE int cmp_br_br(mpr*, mpr*); + +/* + * Initialization + */ +static n_real zero, one, two; + +static char *fatal_error_strings[] = { +#define DIVIDE_BY_ZERO 0 + "divide by zero", +#define FLOATING_POINT_OVERFLOW 1 + "floating point overflow", +#define FLOATING_POINT_EXCEPTION 2 + "floating point exception" +}; + +static char *fatal_object_error_strings[] = { +#define NOT_A_NUMBER 0 + "is not a number", +#define NOT_A_REAL_NUMBER 1 + "is not a real number", +#define NOT_AN_INTEGER 2 + "is not an integer" +}; + +/* + * Implementation + */ +static void +fatal_error(int num) +{ + LispDestroy(fatal_error_strings[num]); +} + +static void +fatal_object_error(LispObj *obj, int num) +{ + LispDestroy("%s %s", STROBJ(obj), fatal_object_error_strings[num]); +} + +static void +fatal_builtin_object_error(LispBuiltin *builtin, LispObj *obj, int num) +{ + LispDestroy("%s: %s %s", STRFUN(builtin), STROBJ(obj), + fatal_object_error_strings[num]); +} + +static void +number_init(void) +{ + zero.type = one.type = two.type = N_FIXNUM; + zero.data.fixnum = 0; + one.data.fixnum = 1; + two.data.fixnum = 2; +} + +static double +bi_getd(mpi *bignum) +{ + double value = mpi_getd(bignum); + + if (!finite(value)) + fatal_error(FLOATING_POINT_EXCEPTION); + + return (value); +} + +static double +br_getd(mpr *bigratio) +{ + double value = mpr_getd(bigratio); + + if (!finite(value)) + fatal_error(FLOATING_POINT_EXCEPTION); + + return (value); +} + +static LispObj * +number_pi(void) +{ + LispObj *result; +#ifndef M_PI +#define M_PI 3.14159265358979323846 +#endif + result = DFLOAT(M_PI); + + return (result); +} + +static void +set_real_real(n_real *real, n_real *val) +{ + switch (RTYPE(real) = RTYPE(val)) { + case N_FIXNUM: + RFI(real) = RFI(val); + break; + case N_BIGNUM: + RBI(real) = XALLOC(mpi); + mpi_init(RBI(real)); + mpi_set(RBI(real), RBI(val)); + break; + case N_FLONUM: + RFF(real) = RFF(val); + break; + case N_FIXRATIO: + RFRN(real) = RFRN(val); + RFRD(real) = RFRD(val); + break; + case N_BIGRATIO: + RBR(real) = XALLOC(mpr); + mpr_init(RBR(real)); + mpr_set(RBR(real), RBR(val)); + break; + } +} + +static void +set_real_object(n_real *real, LispObj *obj) +{ + switch (OBJECT_TYPE(obj)) { + case LispFixnum_t: + RTYPE(real) = N_FIXNUM; + RFI(real) = OFI(obj); + break; + case LispInteger_t: + RTYPE(real) = N_FIXNUM; + RFI(real) = OII(obj); + break; + case LispBignum_t: + RTYPE(real) = N_BIGNUM; + RBI(real) = XALLOC(mpi); + mpi_init(RBI(real)); + mpi_set(RBI(real), OBI(obj)); + break; + case LispDFloat_t: + RTYPE(real) = N_FLONUM; + RFF(real) = ODF(obj); + break; + case LispRatio_t: + RTYPE(real) = N_FIXRATIO; + RFRN(real) = OFRN(obj); + RFRD(real) = OFRD(obj); + break; + case LispBigratio_t: + RTYPE(real) = N_BIGRATIO; + RBR(real) = XALLOC(mpr); + mpr_init(RBR(real)); + mpr_set(RBR(real), OBR(obj)); + break; + default: + fatal_object_error(obj, NOT_A_REAL_NUMBER); + break; + } +} + +static void +set_number_object(n_number *num, LispObj *obj) +{ + switch (OBJECT_TYPE(obj)) { + case LispFixnum_t: + num->complex = 0; + NRTYPE(num) = N_FIXNUM; + NRFI(num) = OFI(obj); + break; + case LispInteger_t: + num->complex = 0; + NRTYPE(num) = N_FIXNUM; + NRFI(num) = OII(obj); + break; + case LispBignum_t: + num->complex = 0; + NRTYPE(num) = N_BIGNUM; + NRBI(num) = XALLOC(mpi); + mpi_init(NRBI(num)); + mpi_set(NRBI(num), OBI(obj)); + break; + case LispDFloat_t: + num->complex = 0; + NRTYPE(num) = N_FLONUM; + NRFF(num) = ODF(obj); + break; + case LispRatio_t: + num->complex = 0; + NRTYPE(num) = N_FIXRATIO; + NRFRN(num) = OFRN(obj); + NRFRD(num) = OFRD(obj); + break; + case LispBigratio_t: + num->complex = 0; + NRTYPE(num) = N_BIGRATIO; + NRBR(num) = XALLOC(mpr); + mpr_init(NRBR(num)); + mpr_set(NRBR(num), OBR(obj)); + break; + case LispComplex_t: + num->complex = 1; + set_real_object(NREAL(num), OCXR(obj)); + set_real_object(NIMAG(num), OCXI(obj)); + break; + default: + fatal_object_error(obj, NOT_A_NUMBER); + break; + } +} + +static void +clear_real(n_real *real) +{ + if (RTYPE(real) == N_BIGNUM) { + mpi_clear(RBI(real)); + XFREE(RBI(real)); + } + else if (RTYPE(real) == N_BIGRATIO) { + mpr_clear(RBR(real)); + XFREE(RBR(real)); + } +} + +static void +clear_number(n_number *num) +{ + clear_real(NREAL(num)); + if (num->complex) + clear_real(NIMAG(num)); +} + +static LispObj * +make_real_object(n_real *real) +{ + LispObj *obj; + + switch (RTYPE(real)) { + case N_FIXNUM: + if (RFI(real) > MOST_POSITIVE_FIXNUM || + RFI(real) < MOST_NEGATIVE_FIXNUM) { + obj = LispNew(NIL, NIL); + obj->type = LispInteger_t; + OII(obj) = RFI(real); + } + else + obj = FIXNUM(RFI(real)); + break; + case N_BIGNUM: + obj = BIGNUM(RBI(real)); + break; + case N_FLONUM: + obj = DFLOAT(RFF(real)); + break; + case N_FIXRATIO: + obj = LispNew(NIL, NIL); + obj->type = LispRatio_t; + OFRN(obj) = RFRN(real); + OFRD(obj) = RFRD(real); + break; + case N_BIGRATIO: + obj = BIGRATIO(RBR(real)); + break; + default: + obj = NIL; + break; + } + + return (obj); +} + +static LispObj * +make_number_object(n_number *num) +{ + LispObj *obj; + + if (num->complex) { + GC_ENTER(); + + obj = LispNew(NIL, NIL); + GC_PROTECT(obj); + OCXI(obj) = NIL; + obj->type = LispComplex_t; + OCXR(obj) = make_real_object(NREAL(num)); + OCXI(obj) = make_real_object(NIMAG(num)); + GC_LEAVE(); + } + else { + switch (NRTYPE(num)) { + case N_FIXNUM: + if (NRFI(num) > MOST_POSITIVE_FIXNUM || + NRFI(num) < MOST_NEGATIVE_FIXNUM) { + obj = LispNew(NIL, NIL); + obj->type = LispInteger_t; + OII(obj) = NRFI(num); + } + else + obj = FIXNUM(NRFI(num)); + break; + case N_BIGNUM: + obj = BIGNUM(NRBI(num)); + break; + case N_FLONUM: + obj = DFLOAT(NRFF(num)); + break; + case N_FIXRATIO: + obj = LispNew(NIL, NIL); + obj->type = LispRatio_t; + OFRN(obj) = NRFRN(num); + OFRD(obj) = NRFRD(num); + break; + case N_BIGRATIO: + obj = BIGRATIO(NRBR(num)); + break; + default: + obj = NIL; + break; + } + } + + return (obj); +} + +#define DEFOP_REAL_REAL(OP) \ +OP##_real_real(n_real *real, n_real *val) \ +{ \ + switch (RTYPE(real)) { \ + case N_FIXNUM: \ + switch (RTYPE(val)) { \ + case N_FIXNUM: \ + r##OP##_fi_fi(real, RFI(val)); \ + break; \ + case N_BIGNUM: \ + r##OP##_fi_bi(real, RBI(val)); \ + break; \ + case N_FLONUM: \ + r##OP##_flonum(real, (double)RFI(real), RFF(val)); \ + break; \ + case N_FIXRATIO: \ + r##OP##_fi_fr(real, RFRN(val), RFRD(val)); \ + break; \ + case N_BIGRATIO: \ + r##OP##_fi_br(real, RBR(val)); \ + break; \ + } \ + break; \ + case N_BIGNUM: \ + switch (RTYPE(val)) { \ + case N_FIXNUM: \ + r##OP##_bi_fi(real, RFI(val)); \ + break; \ + case N_BIGNUM: \ + r##OP##_bi_bi(real, RBI(val)); \ + break; \ + case N_FLONUM: \ + r##OP##_flonum(real, bi_getd(RBI(real)), RFF(val)); \ + break; \ + case N_FIXRATIO: \ + r##OP##_bi_fr(real, RFRN(val), RFRD(val)); \ + break; \ + case N_BIGRATIO: \ + r##OP##_bi_br(real, RBR(val)); \ + break; \ + } \ + break; \ + case N_FLONUM: \ + switch (RTYPE(val)) { \ + case N_FIXNUM: \ + r##OP##_flonum(real, RFF(real), (double)RFI(val)); \ + break; \ + case N_BIGNUM: \ + r##OP##_flonum(real, RFF(real), bi_getd(RBI(val))); \ + break; \ + case N_FLONUM: \ + r##OP##_flonum(real, RFF(real), RFF(val)); \ + break; \ + case N_FIXRATIO: \ + r##OP##_flonum(real, RFF(real), \ + (double)RFRN(val) / (double)RFRD(val));\ + break; \ + case N_BIGRATIO: \ + r##OP##_flonum(real, RFF(real), br_getd(RBR(val))); \ + break; \ + } \ + break; \ + case N_FIXRATIO: \ + switch (RTYPE(val)) { \ + case N_FIXNUM: \ + r##OP##_fr_fi(real, RFI(val)); \ + break; \ + case N_BIGNUM: \ + r##OP##_fr_bi(real, RBI(val)); \ + break; \ + case N_FLONUM: \ + r##OP##_flonum(real, \ + (double)RFRN(real) / (double)RFRD(real),\ + RFF(val)); \ + break; \ + case N_FIXRATIO: \ + r##OP##_fr_fr(real, RFRN(val), RFRD(val)); \ + break; \ + case N_BIGRATIO: \ + r##OP##_fr_br(real, RBR(val)); \ + break; \ + } \ + break; \ + case N_BIGRATIO: \ + switch (RTYPE(val)) { \ + case N_FIXNUM: \ + r##OP##_br_fi(real, RFI(val)); \ + break; \ + case N_BIGNUM: \ + r##OP##_br_bi(real, RBI(val)); \ + break; \ + case N_FLONUM: \ + r##OP##_flonum(real, br_getd(RBR(real)), RFF(val)); \ + break; \ + case N_FIXRATIO: \ + r##OP##_br_fr(real, RFRN(val), RFRD(val)); \ + break; \ + case N_BIGRATIO: \ + r##OP##_br_br(real, RBR(val)); \ + break; \ + } \ + break; \ + } \ +} + +static void +DEFOP_REAL_REAL(add) + +static void +DEFOP_REAL_REAL(sub) + +static void +DEFOP_REAL_REAL(div) + +static void +DEFOP_REAL_REAL(mul) + + +#define DEFOP_REAL_OBJECT(OP) \ +OP##_real_object(n_real *real, LispObj *obj) \ +{ \ + switch (OBJECT_TYPE(obj)) { \ + case LispFixnum_t: \ + switch (RTYPE(real)) { \ + case N_FIXNUM: \ + r##OP##_fi_fi(real, OFI(obj)); \ + break; \ + case N_BIGNUM: \ + r##OP##_bi_fi(real, OFI(obj)); \ + break; \ + case N_FLONUM: \ + r##OP##_flonum(real, RFF(real), (double)OFI(obj)); \ + break; \ + case N_FIXRATIO: \ + r##OP##_fr_fi(real, OFI(obj)); \ + break; \ + case N_BIGRATIO: \ + r##OP##_br_fi(real, OFI(obj)); \ + break; \ + } \ + break; \ + case LispInteger_t: \ + switch (RTYPE(real)) { \ + case N_FIXNUM: \ + r##OP##_fi_fi(real, OII(obj)); \ + break; \ + case N_BIGNUM: \ + r##OP##_bi_fi(real, OII(obj)); \ + break; \ + case N_FLONUM: \ + r##OP##_flonum(real, RFF(real), (double)OII(obj)); \ + break; \ + case N_FIXRATIO: \ + r##OP##_fr_fi(real, OII(obj)); \ + break; \ + case N_BIGRATIO: \ + r##OP##_br_fi(real, OII(obj)); \ + break; \ + } \ + break; \ + case LispBignum_t: \ + switch (RTYPE(real)) { \ + case N_FIXNUM: \ + r##OP##_fi_bi(real, OBI(obj)); \ + break; \ + case N_BIGNUM: \ + r##OP##_bi_bi(real, OBI(obj)); \ + break; \ + case N_FLONUM: \ + r##OP##_flonum(real, RFF(real), bi_getd(OBI(obj))); \ + break; \ + case N_FIXRATIO: \ + r##OP##_fr_bi(real, OBI(obj)); \ + break; \ + case N_BIGRATIO: \ + r##OP##_br_bi(real, OBI(obj)); \ + break; \ + } \ + break; \ + case LispDFloat_t: \ + switch (RTYPE(real)) { \ + case N_FIXNUM: \ + r##OP##_flonum(real, (double)RFI(real), ODF(obj)); \ + break; \ + case N_BIGNUM: \ + r##OP##_flonum(real, bi_getd(RBI(real)), ODF(obj)); \ + break; \ + case N_FLONUM: \ + r##OP##_flonum(real, RFF(real), ODF(obj)); \ + break; \ + case N_FIXRATIO: \ + r##OP##_flonum(real, \ + (double)RFRN(real) / (double)RFRD(real),\ + ODF(obj)); \ + break; \ + case N_BIGRATIO: \ + r##OP##_flonum(real, br_getd(RBR(real)), ODF(obj)); \ + break; \ + } \ + break; \ + case LispRatio_t: \ + switch (RTYPE(real)) { \ + case N_FIXNUM: \ + r##OP##_fi_fr(real, OFRN(obj), OFRD(obj)); \ + break; \ + case N_BIGNUM: \ + r##OP##_bi_fr(real, OFRN(obj), OFRD(obj)); \ + break; \ + case N_FLONUM: \ + r##OP##_flonum(real, RFF(real), \ + (double)OFRN(obj) / (double)OFRD(obj)); \ + break; \ + case N_FIXRATIO: \ + r##OP##_fr_fr(real, OFRN(obj), OFRD(obj)); \ + break; \ + case N_BIGRATIO: \ + r##OP##_br_fr(real, OFRN(obj), OFRD(obj)); \ + break; \ + } \ + break; \ + case LispBigratio_t: \ + switch (RTYPE(real)) { \ + case N_FIXNUM: \ + r##OP##_fi_br(real, OBR(obj)); \ + break; \ + case N_BIGNUM: \ + r##OP##_bi_br(real, OBR(obj)); \ + break; \ + case N_FLONUM: \ + r##OP##_flonum(real, RFF(real), br_getd(OBR(obj))); \ + break; \ + case N_FIXRATIO: \ + r##OP##_fr_br(real, OBR(obj)); \ + break; \ + case N_BIGRATIO: \ + r##OP##_br_br(real, OBR(obj)); \ + break; \ + } \ + break; \ + default: \ + fatal_object_error(obj, NOT_A_REAL_NUMBER); \ + break; \ + } \ +} + +static void +DEFOP_REAL_OBJECT(add) + +static void +DEFOP_REAL_OBJECT(sub) + +static void +DEFOP_REAL_OBJECT(div) + +static void +DEFOP_REAL_OBJECT(mul) + + +#define DEFOP_NUMBER_OBJECT(OP) \ +OP##_number_object(n_number *num, LispObj *obj) \ +{ \ + if (num->complex) { \ + switch (OBJECT_TYPE(obj)) { \ + case LispFixnum_t: \ + case LispInteger_t: \ + case LispBignum_t: \ + case LispDFloat_t: \ + case LispRatio_t: \ + case LispBigratio_t: \ + n##OP##_cx_re(num, obj); \ + break; \ + case LispComplex_t: \ + n##OP##_cx_cx(num, obj); \ + break; \ + default: \ + fatal_object_error(obj, NOT_A_NUMBER); \ + break; \ + } \ + } \ + else { \ + switch (OBJECT_TYPE(obj)) { \ + case LispFixnum_t: \ + switch (NRTYPE(num)) { \ + case N_FIXNUM: \ + r##OP##_fi_fi(NREAL(num), OFI(obj)); \ + break; \ + case N_BIGNUM: \ + r##OP##_bi_fi(NREAL(num), OFI(obj)); \ + break; \ + case N_FLONUM: \ + r##OP##_flonum(NREAL(num), NRFF(num), \ + (double)OFI(obj)); \ + break; \ + case N_FIXRATIO: \ + r##OP##_fr_fi(NREAL(num), OFI(obj)); \ + break; \ + case N_BIGRATIO: \ + r##OP##_br_fi(NREAL(num), OFI(obj)); \ + break; \ + } \ + break; \ + case LispInteger_t: \ + switch (NRTYPE(num)) { \ + case N_FIXNUM: \ + r##OP##_fi_fi(NREAL(num), OII(obj)); \ + break; \ + case N_BIGNUM: \ + r##OP##_bi_fi(NREAL(num), OII(obj)); \ + break; \ + case N_FLONUM: \ + r##OP##_flonum(NREAL(num), NRFF(num), \ + (double)OII(obj)); \ + break; \ + case N_FIXRATIO: \ + r##OP##_fr_fi(NREAL(num), OII(obj)); \ + break; \ + case N_BIGRATIO: \ + r##OP##_br_fi(NREAL(num), OII(obj)); \ + break; \ + } \ + break; \ + case LispBignum_t: \ + switch (NRTYPE(num)) { \ + case N_FIXNUM: \ + r##OP##_fi_bi(NREAL(num), OBI(obj)); \ + break; \ + case N_BIGNUM: \ + r##OP##_bi_bi(NREAL(num), OBI(obj)); \ + break; \ + case N_FLONUM: \ + r##OP##_flonum(NREAL(num), NRFF(num), \ + bi_getd(OBI(obj))); \ + break; \ + case N_FIXRATIO: \ + r##OP##_fr_bi(NREAL(num), OBI(obj)); \ + break; \ + case N_BIGRATIO: \ + r##OP##_br_bi(NREAL(num), OBI(obj)); \ + break; \ + } \ + break; \ + case LispDFloat_t: \ + switch (NRTYPE(num)) { \ + case N_FIXNUM: \ + r##OP##_flonum(NREAL(num), (double)NRFI(num), \ + ODF(obj)); \ + break; \ + case N_BIGNUM: \ + r##OP##_flonum(NREAL(num), bi_getd(NRBI(num)), \ + ODF(obj)); \ + break; \ + case N_FLONUM: \ + r##OP##_flonum(NREAL(num), NRFF(num), ODF(obj));\ + break; \ + case N_FIXRATIO: \ + r##OP##_flonum(NREAL(num), \ + (double)NRFRN(num) / \ + (double)NRFRD(num), \ + ODF(obj)); \ + break; \ + case N_BIGRATIO: \ + r##OP##_flonum(NREAL(num), br_getd(NRBR(num)), \ + ODF(obj)); \ + break; \ + } \ + break; \ + case LispRatio_t: \ + switch (NRTYPE(num)) { \ + case N_FIXNUM: \ + r##OP##_fi_fr(NREAL(num), OFRN(obj), OFRD(obj));\ + break; \ + case N_BIGNUM: \ + r##OP##_bi_fr(NREAL(num), OFRN(obj), OFRD(obj));\ + break; \ + case N_FLONUM: \ + r##OP##_flonum(NREAL(num), NRFF(num), \ + (double)OFRN(obj) / \ + (double)OFRD(obj)); \ + break; \ + case N_FIXRATIO: \ + r##OP##_fr_fr(NREAL(num), OFRN(obj), OFRD(obj));\ + break; \ + case N_BIGRATIO: \ + r##OP##_br_fr(NREAL(num), OFRN(obj), OFRD(obj));\ + break; \ + } \ + break; \ + case LispBigratio_t: \ + switch (NRTYPE(num)) { \ + case N_FIXNUM: \ + r##OP##_fi_br(NREAL(num), OBR(obj)); \ + break; \ + case N_BIGNUM: \ + r##OP##_bi_br(NREAL(num), OBR(obj)); \ + break; \ + case N_FLONUM: \ + r##OP##_flonum(NREAL(num), NRFF(num), \ + br_getd(OBR(obj))); \ + break; \ + case N_FIXRATIO: \ + r##OP##_fr_br(NREAL(num), OBR(obj)); \ + break; \ + case N_BIGRATIO: \ + r##OP##_br_br(NREAL(num), OBR(obj)); \ + break; \ + } \ + break; \ + case LispComplex_t: \ + n##OP##_re_cx(num, obj); \ + break; \ + default: \ + fatal_object_error(obj, NOT_A_NUMBER); \ + break; \ + } \ + } \ +} + +static void +DEFOP_NUMBER_OBJECT(add) + +static void +DEFOP_NUMBER_OBJECT(sub) + +static void +DEFOP_NUMBER_OBJECT(div) + +static void +DEFOP_NUMBER_OBJECT(mul) + + +/************************************************************************ + * ABS + ************************************************************************/ +static void +abs_real(n_real *real) +{ + switch (RTYPE(real)) { + case N_FIXNUM: rabs_fi(real); break; + case N_BIGNUM: rabs_bi(real); break; + case N_FLONUM: rabs_ff(real); break; + case N_FIXRATIO: rabs_fr(real); break; + case N_BIGRATIO: rabs_br(real); break; + } +} + +static void +abs_number(n_number *num) +{ + if (num->complex) + nabs_cx(num); + else { + switch (NRTYPE(num)) { + case N_FIXNUM: rabs_fi(NREAL(num)); break; + case N_BIGNUM: rabs_bi(NREAL(num)); break; + case N_FLONUM: rabs_ff(NREAL(num)); break; + case N_FIXRATIO: rabs_fr(NREAL(num)); break; + case N_BIGRATIO: rabs_br(NREAL(num)); break; + } + } +} + +static void +nabs_cx(n_number *num) +{ + n_real temp; + + abs_real(NREAL(num)); + abs_real(NIMAG(num)); + + if (cmp_real_real(NREAL(num), NIMAG(num)) < 0) { + memcpy(&temp, NIMAG(num), sizeof(n_real)); + memcpy(NIMAG(num), NREAL(num), sizeof(n_real)); + memcpy(NREAL(num), &temp, sizeof(n_real)); + } + + if (cmp_real_real(NIMAG(num), &zero) == 0) { + num->complex = 0; + if (NITYPE(num) == N_FLONUM) { + /* change number type */ + temp.type = N_FLONUM; + temp.data.flonum = 1.0; + mul_real_real(NREAL(num), &temp); + } + else + clear_real(NIMAG(num)); + } + else { + div_real_real(NIMAG(num), NREAL(num)); + set_real_real(&temp, NIMAG(num)); + mul_real_real(NIMAG(num), &temp); + clear_real(&temp); + + add_real_real(NIMAG(num), &one); + sqrt_real(NIMAG(num)); + + mul_real_real(NIMAG(num), NREAL(num)); + clear_real(NREAL(num)); + memcpy(NREAL(num), NIMAG(num), sizeof(n_real)); + num->complex = 0; + } +} + +static INLINE void +rabs_fi(n_real *real) +{ + if (RFI(real) < 0) + rneg_fi(real); +} + +static INLINE void +rabs_bi(n_real *real) +{ + if (mpi_cmpi(RBI(real), 0) < 0) + mpi_neg(RBI(real), RBI(real)); +} + +static INLINE void +rabs_ff(n_real *real) +{ + if (RFF(real) < 0.0) + RFF(real) = -RFF(real); +} + +static INLINE void +rabs_fr(n_real *real) +{ + if (RFRN(real) < 0) + rneg_fr(real); +} + +static INLINE void +rabs_br(n_real *real) +{ + if (mpi_cmpi(RBRN(real), 0) < 0) + mpi_neg(RBRN(real), RBRN(real)); +} + + +/************************************************************************ + * NEG + ************************************************************************/ +static void +neg_real(n_real *real) +{ + switch (RTYPE(real)) { + case N_FIXNUM: rneg_fi(real); break; + case N_BIGNUM: rneg_bi(real); break; + case N_FLONUM: rneg_ff(real); break; + case N_FIXRATIO: rneg_fr(real); break; + case N_BIGRATIO: rneg_br(real); break; + } +} + +static void +neg_number(n_number *num) +{ + if (num->complex) { + neg_real(NREAL(num)); + neg_real(NIMAG(num)); + } + else { + switch (NRTYPE(num)) { + case N_FIXNUM: rneg_fi(NREAL(num)); break; + case N_BIGNUM: rneg_bi(NREAL(num)); break; + case N_FLONUM: rneg_ff(NREAL(num)); break; + case N_FIXRATIO: rneg_fr(NREAL(num)); break; + case N_BIGRATIO: rneg_br(NREAL(num)); break; + } + } +} + +static void +rneg_fi(n_real *real) +{ + if (RFI(real) == MINSLONG) { + mpi *bigi = XALLOC(mpi); + + mpi_init(bigi); + mpi_seti(bigi, RFI(real)); + mpi_neg(bigi, bigi); + RTYPE(real) = N_BIGNUM; + RBI(real) = bigi; + } + else + RFI(real) = -RFI(real); +} + +static INLINE void +rneg_bi(n_real *real) +{ + mpi_neg(RBI(real), RBI(real)); +} + +static INLINE void +rneg_ff(n_real *real) +{ + RFF(real) = -RFF(real); +} + +static void +rneg_fr(n_real *real) +{ + if (RFRN(real) == MINSLONG) { + mpr *bigr = XALLOC(mpr); + + mpr_init(bigr); + mpr_seti(bigr, RFRN(real), RFRD(real)); + mpi_neg(mpr_num(bigr), mpr_num(bigr)); + RTYPE(real) = N_BIGRATIO; + RBR(real) = bigr; + } + else + RFRN(real) = -RFRN(real); +} + +static INLINE void +rneg_br(n_real *real) +{ + mpi_neg(RBRN(real), RBRN(real)); +} + + +/************************************************************************ + * SQRT + ************************************************************************/ +static void +sqrt_real(n_real *real) +{ + switch (RTYPE(real)) { + case N_FIXNUM: + case N_BIGNUM: + rsqrt_xi(real); + break; + case N_FLONUM: + rsqrt_ff(real); + break; + case N_FIXRATIO: + case N_BIGRATIO: + rsqrt_xr(real); + break; + } +} + +static void +sqrt_number(n_number *num) +{ + if (num->complex) + nsqrt_cx(num); + else { + switch (NRTYPE(num)) { + case N_FIXNUM: + case N_BIGNUM: + nsqrt_xi(num); + break; + case N_FLONUM: + nsqrt_ff(num); + break; + case N_FIXRATIO: + case N_BIGRATIO: + nsqrt_xr(num); + break; + } + } +} + +static void +rsqrt_xi(n_real *real) +{ + int exact; + mpi bignum; + + if (cmp_real_real(real, &zero) < 0) + fatal_error(FLOATING_POINT_EXCEPTION); + + mpi_init(&bignum); + if (RTYPE(real) == N_BIGNUM) + exact = mpi_sqrt(&bignum, RBI(real)); + else { + mpi tmp; + + mpi_init(&tmp); + mpi_seti(&tmp, RFI(real)); + exact = mpi_sqrt(&bignum, &tmp); + mpi_clear(&tmp); + } + if (exact) { + if (RTYPE(real) == N_BIGNUM) { + mpi_set(RBI(real), &bignum); + rbi_canonicalize(real); + } + else + RFI(real) = mpi_geti(&bignum); + } + else { + double value; + + if (RTYPE(real) == N_BIGNUM) { + value = bi_getd(RBI(real)); + RCLEAR_BI(real); + } + else + value = (double)RFI(real); + + value = sqrt(value); + RTYPE(real) = N_FLONUM; + RFF(real) = value; + } + mpi_clear(&bignum); +} + +static void +rsqrt_xr(n_real *real) +{ + n_real num, den; + + if (cmp_real_real(real, &zero) < 0) + fatal_error(FLOATING_POINT_EXCEPTION); + + if (RTYPE(real) == N_FIXRATIO) { + num.type = den.type = N_FIXNUM; + num.data.fixnum = RFRN(real); + den.data.fixnum = RFRD(real); + } + else { + mpi *bignum; + + if (mpi_fiti(RBRN(real))) { + num.type = N_FIXNUM; + num.data.fixnum = mpi_geti(RBRN(real)); + } + else { + bignum = XALLOC(mpi); + mpi_init(bignum); + mpi_set(bignum, RBRN(real)); + num.type = N_BIGNUM; + num.data.bignum = bignum; + } + + if (mpi_fiti(RBRD(real))) { + den.type = N_FIXNUM; + den.data.fixnum = mpi_geti(RBRD(real)); + } + else { + bignum = XALLOC(mpi); + mpi_init(bignum); + mpi_set(bignum, RBRD(real)); + den.type = N_BIGNUM; + den.data.bignum = bignum; + } + } + + rsqrt_xi(&num); + rsqrt_xi(&den); + + clear_real(real); + memcpy(real, &num, sizeof(n_real)); + div_real_real(real, &den); + clear_real(&den); +} + +static void +rsqrt_ff(n_real *real) +{ + if (RFF(real) < 0.0) + fatal_error(FLOATING_POINT_EXCEPTION); + RFF(real) = sqrt(RFF(real)); +} + + +static void +nsqrt_cx(n_number *num) +{ + n_number mag; + n_real *real, *imag; + + real = &(mag.real); + imag = &(mag.imag); + set_real_real(real, NREAL(num)); + set_real_real(imag, NIMAG(num)); + mag.complex = 1; + + nabs_cx(&mag); /* this will free the imag part data */ + if (cmp_real_real(real, &zero) == 0) { + clear_number(num); + memcpy(NREAL(num), real, sizeof(n_real)); + clear_real(real); + num->complex = 0; + return; + } + else if (cmp_real_real(NREAL(num), &zero) > 0) { + /* R = sqrt((mag + Ra) / 2) */ + add_real_real(NREAL(num), real); + clear_real(real); + div_real_real(NREAL(num), &two); + sqrt_real(NREAL(num)); + + /* I = Ia / R / 2 */ + div_real_real(NIMAG(num), NREAL(num)); + div_real_real(NIMAG(num), &two); + } + else { + /* remember old imag part */ + memcpy(imag, NIMAG(num), sizeof(n_real)); + + /* I = sqrt((mag - Ra) / 2) */ + memcpy(NIMAG(num), real, sizeof(n_real)); + sub_real_real(NIMAG(num), NREAL(num)); + div_real_real(NIMAG(num), &two); + sqrt_real(NIMAG(num)); + if (cmp_real_real(imag, &zero) < 0) + neg_real(NIMAG(num)); + + /* R = Ia / I / 2 */ + clear_real(NREAL(num)); + /* start with old imag part */ + memcpy(NREAL(num), imag, sizeof(n_real)); + div_real_real(NREAL(num), NIMAG(num)); + div_real_real(NREAL(num), &two); + } + + ncx_canonicalize(num); +} + +static void +nsqrt_xi(n_number *num) +{ + if (cmp_real_real(NREAL(num), &zero) < 0) { + memcpy(NIMAG(num), NREAL(num), sizeof(n_real)); + neg_real(NIMAG(num)); + rsqrt_xi(NIMAG(num)); + NRTYPE(num) = N_FIXNUM; + NRFI(num) = 0; + num->complex = 1; + } + else + rsqrt_xi(NREAL(num)); +} + +static void +nsqrt_ff(n_number *num) +{ + double value; + + if (NRFF(num) < 0.0) { + value = sqrt(-NRFF(num)); + + NITYPE(num) = N_FLONUM; + NIFF(num) = value; + NRTYPE(num) = N_FIXNUM; + NRFI(num) = 0; + num->complex = 1; + } + else { + value = sqrt(NRFF(num)); + NRFF(num) = value; + } +} + +static void +nsqrt_xr(n_number *num) +{ + if (cmp_real_real(NREAL(num), &zero) < 0) { + memcpy(NIMAG(num), NREAL(num), sizeof(n_real)); + neg_real(NIMAG(num)); + rsqrt_xr(NIMAG(num)); + NRTYPE(num) = N_FIXNUM; + NRFI(num) = 0; + num->complex = 1; + } + else + rsqrt_xr(NREAL(num)); +} + + +/************************************************************************ + * MOD + ************************************************************************/ +static void +mod_real_real(n_real *real, n_real *val) +{ + /* Assume both operands are integers */ + switch (RTYPE(real)) { + case N_FIXNUM: + switch (RTYPE(val)) { + case N_FIXNUM: + rmod_fi_fi(real, RFI(val)); + break; + case N_BIGNUM: + rmod_fi_bi(real, RBI(val)); + break; + } + break; + case N_BIGNUM: + switch (RTYPE(val)) { + case N_FIXNUM: + rmod_bi_fi(real, RFI(val)); + break; + case N_BIGNUM: + rmod_bi_bi(real, RBI(val)); + break; + } + break; + } +} + +static void +mod_real_object(n_real *real, LispObj *obj) +{ + switch (RTYPE(real)) { + case N_FIXNUM: + switch (OBJECT_TYPE(obj)) { + case LispFixnum_t: + rmod_fi_fi(real, OFI(obj)); + return; + case LispInteger_t: + rmod_fi_fi(real, OII(obj)); + return; + case LispBignum_t: + rmod_fi_bi(real, OBI(obj)); + return; + default: + break; + } + break; + case N_BIGNUM: + switch (OBJECT_TYPE(obj)) { + case LispFixnum_t: + rmod_bi_fi(real, OFI(obj)); + return; + case LispInteger_t: + rmod_bi_fi(real, OII(obj)); + return; + case LispBignum_t: + rmod_bi_bi(real, OBI(obj)); + return; + default: + break; + } + break; + /* Assume the n_real object is an integer */ + } + fatal_object_error(obj, NOT_AN_INTEGER); +} + +static void +rmod_fi_fi(n_real *real, long fi) +{ + if (fi == 0) + fatal_error(DIVIDE_BY_ZERO); + + if ((RFI(real) < 0) ^ (fi < 0)) + RFI(real) = (RFI(real) % fi) + fi; + else if (RFI(real) == MINSLONG || fi == MINSLONG) { + mpi bignum; + + mpi_init(&bignum); + mpi_seti(&bignum, RFI(real)); + RFI(real) = mpi_modi(&bignum, fi); + mpi_clear(&bignum); + } + else + RFI(real) = RFI(real) % fi; +} + +static void +rmod_fi_bi(n_real *real, mpi *bignum) +{ + mpi *bigi; + + if (mpi_cmpi(bignum, 0) == 0) + fatal_error(DIVIDE_BY_ZERO); + + bigi = XALLOC(mpi); + mpi_init(bigi); + mpi_seti(bigi, RFI(real)); + mpi_mod(bigi, bigi, bignum); + RTYPE(real) = N_BIGNUM; + RBI(real) = bigi; + rbi_canonicalize(real); +} + +static void +rmod_bi_fi(n_real *real, long fi) +{ + mpi iop; + + if (fi == 0) + fatal_error(DIVIDE_BY_ZERO); + + mpi_init(&iop); + mpi_seti(&iop, fi); + mpi_mod(RBI(real), RBI(real), &iop); + mpi_clear(&iop); + rbi_canonicalize(real); +} + +static void +rmod_bi_bi(n_real *real, mpi *bignum) +{ + if (mpi_cmpi(bignum, 0) == 0) + fatal_error(DIVIDE_BY_ZERO); + + mpi_mod(RBI(real), RBI(real), bignum); + rbi_canonicalize(real); +} + +/************************************************************************ + * REM + ************************************************************************/ +static void +rem_real_object(n_real *real, LispObj *obj) +{ + switch (RTYPE(real)) { + case N_FIXNUM: + switch (OBJECT_TYPE(obj)) { + case LispFixnum_t: + rrem_fi_fi(real, OFI(obj)); + return; + case LispInteger_t: + rrem_fi_fi(real, OII(obj)); + return; + case LispBignum_t: + rrem_fi_bi(real, OBI(obj)); + return; + default: + break; + } + break; + case N_BIGNUM: + switch (OBJECT_TYPE(obj)) { + case LispFixnum_t: + rrem_bi_fi(real, OFI(obj)); + return; + case LispInteger_t: + rrem_bi_fi(real, OII(obj)); + return; + case LispBignum_t: + rrem_bi_bi(real, OBI(obj)); + return; + default: + break; + } + break; + /* Assume the n_real object is an integer */ + } + fatal_object_error(obj, NOT_AN_INTEGER); +} + +static void +rrem_fi_fi(n_real *real, long fi) +{ + if (fi == 0) + fatal_error(DIVIDE_BY_ZERO); + + if (RFI(real) == MINSLONG || fi == MINSLONG) { + mpi bignum; + + mpi_init(&bignum); + mpi_seti(&bignum, RFI(real)); + RFI(real) = mpi_remi(&bignum, fi); + mpi_clear(&bignum); + } + else + RFI(real) = RFI(real) % fi; +} + +static void +rrem_fi_bi(n_real *real, mpi *bignum) +{ + mpi *bigi; + + if (mpi_cmpi(bignum, 0) == 0) + fatal_error(DIVIDE_BY_ZERO); + + bigi = XALLOC(mpi); + mpi_init(bigi); + mpi_seti(bigi, RFI(real)); + mpi_rem(bigi, bigi, bignum); + RTYPE(real) = N_BIGNUM; + RBI(real) = bigi; + rbi_canonicalize(real); +} + +static void +rrem_bi_fi(n_real *real, long fi) +{ + mpi iop; + + if (fi == 0) + fatal_error(DIVIDE_BY_ZERO); + + mpi_init(&iop); + mpi_seti(&iop, fi); + mpi_rem(RBI(real), RBI(real), &iop); + mpi_clear(&iop); + rbi_canonicalize(real); +} + +static void +rrem_bi_bi(n_real *real, mpi *bignum) +{ + if (mpi_cmpi(bignum, 0) == 0) + fatal_error(DIVIDE_BY_ZERO); + + mpi_rem(RBI(real), RBI(real), bignum); + rbi_canonicalize(real); +} + + +/************************************************************************ + * GCD + ************************************************************************/ +static void +gcd_real_object(n_real *real, LispObj *obj) +{ + if (!INTEGERP(obj)) + fatal_object_error(obj, NOT_AN_INTEGER); + + /* check for zero operand */ + if (cmp_real_real(real, &zero) == 0) + set_real_object(real, obj); + else if (cmp_real_object(&zero, obj) != 0) { + n_real rest, temp; + + set_real_object(&rest, obj); + for (;;) { + mod_real_real(&rest, real); + if (cmp_real_real(&rest, &zero) == 0) + break; + memcpy(&temp, real, sizeof(n_real)); + memcpy(real, &rest, sizeof(n_real)); + memcpy(&rest, &temp, sizeof(n_real)); + } + clear_real(&rest); + } +} + +/************************************************************************ + * AND + ************************************************************************/ +static void +and_real_object(n_real *real, LispObj *obj) +{ + mpi *bigi, iop; + + switch (OBJECT_TYPE(obj)) { + case LispFixnum_t: + switch (RTYPE(real)) { + case N_FIXNUM: + RFI(real) &= OFI(obj); + break; + case N_BIGNUM: + mpi_init(&iop); + mpi_seti(&iop, OFI(obj)); + mpi_and(RBI(real), RBI(real), &iop); + mpi_clear(&iop); + rbi_canonicalize(real); + break; + } + break; + case LispInteger_t: + switch (RTYPE(real)) { + case N_FIXNUM: + RFI(real) &= OII(obj); + break; + case N_BIGNUM: + mpi_init(&iop); + mpi_seti(&iop, OII(obj)); + mpi_and(RBI(real), RBI(real), &iop); + mpi_clear(&iop); + rbi_canonicalize(real); + break; + } + break; + case LispBignum_t: + switch (RTYPE(real)) { + case N_FIXNUM: + bigi = XALLOC(mpi); + mpi_init(bigi); + mpi_seti(bigi, RFI(real)); + mpi_and(bigi, bigi, OBI(obj)); + RTYPE(real) = N_BIGNUM; + RBI(real) = bigi; + rbi_canonicalize(real); + break; + case N_BIGNUM: + mpi_and(RBI(real), RBI(real), OBI(obj)); + rbi_canonicalize(real); + break; + } + break; + default: + fatal_object_error(obj, NOT_AN_INTEGER); + break; + } +} + + +/************************************************************************ + * EQV + ************************************************************************/ +static void +eqv_real_object(n_real *real, LispObj *obj) +{ + mpi *bigi, iop; + + switch (OBJECT_TYPE(obj)) { + case LispFixnum_t: + switch (RTYPE(real)) { + case N_FIXNUM: + RFI(real) ^= ~OFI(obj); + break; + case N_BIGNUM: + mpi_init(&iop); + mpi_seti(&iop, OFI(obj)); + mpi_com(&iop, &iop); + mpi_xor(RBI(real), RBI(real), &iop); + mpi_clear(&iop); + rbi_canonicalize(real); + break; + } + break; + case LispInteger_t: + switch (RTYPE(real)) { + case N_FIXNUM: + RFI(real) ^= ~OII(obj); + break; + case N_BIGNUM: + mpi_init(&iop); + mpi_seti(&iop, OII(obj)); + mpi_com(&iop, &iop); + mpi_xor(RBI(real), RBI(real), &iop); + mpi_clear(&iop); + rbi_canonicalize(real); + break; + } + break; + case LispBignum_t: + switch (RTYPE(real)) { + case N_FIXNUM: + bigi = XALLOC(mpi); + mpi_init(bigi); + mpi_seti(bigi, RFI(real)); + mpi_com(bigi, bigi); + mpi_xor(bigi, bigi, OBI(obj)); + RTYPE(real) = N_BIGNUM; + RBI(real) = bigi; + rbi_canonicalize(real); + break; + case N_BIGNUM: + mpi_com(RBI(real), RBI(real)); + mpi_xor(RBI(real), RBI(real), OBI(obj)); + rbi_canonicalize(real); + break; + } + break; + default: + fatal_object_error(obj, NOT_AN_INTEGER); + break; + } +} + + +/************************************************************************ + * IOR + ************************************************************************/ +static void +ior_real_object(n_real *real, LispObj *obj) +{ + mpi *bigi, iop; + + switch (OBJECT_TYPE(obj)) { + case LispFixnum_t: + switch (RTYPE(real)) { + case N_FIXNUM: + RFI(real) |= OFI(obj); + break; + case N_BIGNUM: + mpi_init(&iop); + mpi_seti(&iop, OFI(obj)); + mpi_ior(RBI(real), RBI(real), &iop); + mpi_clear(&iop); + rbi_canonicalize(real); + break; + } + break; + case LispInteger_t: + switch (RTYPE(real)) { + case N_FIXNUM: + RFI(real) |= OII(obj); + break; + case N_BIGNUM: + mpi_init(&iop); + mpi_seti(&iop, OII(obj)); + mpi_ior(RBI(real), RBI(real), &iop); + mpi_clear(&iop); + rbi_canonicalize(real); + break; + } + break; + case LispBignum_t: + switch (RTYPE(real)) { + case N_FIXNUM: + bigi = XALLOC(mpi); + mpi_init(bigi); + mpi_seti(bigi, RFI(real)); + mpi_ior(bigi, bigi, OBI(obj)); + RTYPE(real) = N_BIGNUM; + RBI(real) = bigi; + rbi_canonicalize(real); + break; + case N_BIGNUM: + mpi_ior(RBI(real), RBI(real), OBI(obj)); + rbi_canonicalize(real); + break; + } + break; + default: + fatal_object_error(obj, NOT_AN_INTEGER); + break; + } +} + + +/************************************************************************ + * NOT + ************************************************************************/ +static void +not_real(n_real *real) +{ + if (RTYPE(real) == N_FIXNUM) + RFI(real) = ~RFI(real); + else { + mpi_com(RBI(real), RBI(real)); + rbi_canonicalize(real); + } +} + +/************************************************************************ + * XOR + ************************************************************************/ +static void +xor_real_object(n_real *real, LispObj *obj) +{ + mpi *bigi, iop; + + switch (OBJECT_TYPE(obj)) { + case LispFixnum_t: + switch (RTYPE(real)) { + case N_FIXNUM: + RFI(real) ^= OFI(obj); + break; + case N_BIGNUM: + mpi_init(&iop); + mpi_seti(&iop, OFI(obj)); + mpi_xor(RBI(real), RBI(real), &iop); + mpi_clear(&iop); + rbi_canonicalize(real); + break; + } + break; + case LispInteger_t: + switch (RTYPE(real)) { + case N_FIXNUM: + RFI(real) ^= OII(obj); + break; + case N_BIGNUM: + mpi_init(&iop); + mpi_seti(&iop, OII(obj)); + mpi_xor(RBI(real), RBI(real), &iop); + mpi_clear(&iop); + rbi_canonicalize(real); + break; + } + break; + case LispBignum_t: + switch (RTYPE(real)) { + case N_FIXNUM: + bigi = XALLOC(mpi); + mpi_init(bigi); + mpi_seti(bigi, RFI(real)); + mpi_xor(bigi, bigi, OBI(obj)); + RTYPE(real) = N_BIGNUM; + RBI(real) = bigi; + rbi_canonicalize(real); + break; + case N_BIGNUM: + mpi_xor(RBI(real), RBI(real), OBI(obj)); + rbi_canonicalize(real); + break; + } + break; + default: + fatal_object_error(obj, NOT_AN_INTEGER); + break; + } +} + + +/************************************************************************ + * DIVIDE + ************************************************************************/ +static void +divide_number_object(n_number *num, LispObj *obj, int fun, int flo) +{ + switch (OBJECT_TYPE(obj)) { + case LispFixnum_t: + switch (NRTYPE(num)) { + case N_FIXNUM: + ndivide_fi_fi(num, OFI(obj), fun, flo); + break; + case N_BIGNUM: + ndivide_xi_xi(num, obj, fun, flo); + break; + case N_FLONUM: + ndivide_flonum(num, NRFF(num), (double)OFI(obj), fun, flo); + break; + case N_FIXRATIO: + case N_BIGRATIO: + ndivide_xr_xi(num, obj, fun, flo); + break; + } + break; + case LispInteger_t: + switch (NRTYPE(num)) { + case N_FIXNUM: + ndivide_fi_fi(num, OII(obj), fun, flo); + break; + case N_BIGNUM: + ndivide_xi_xi(num, obj, fun, flo); + break; + case N_FLONUM: + ndivide_flonum(num, NRFF(num), (double)OII(obj), fun, flo); + break; + case N_FIXRATIO: + case N_BIGRATIO: + ndivide_xr_xi(num, obj, fun, flo); + break; + } + break; + case LispBignum_t: + switch (NRTYPE(num)) { + case N_FIXNUM: + case N_BIGNUM: + ndivide_xi_xi(num, obj, fun, flo); + break; + case N_FLONUM: + ndivide_flonum(num, NRFF(num), bi_getd(OBI(obj)), + fun, flo); + break; + case N_FIXRATIO: + case N_BIGRATIO: + ndivide_xr_xi(num, obj, fun, flo); + break; + } + break; + case LispDFloat_t: + switch (NRTYPE(num)) { + case N_FIXNUM: + ndivide_flonum(num, (double)NRFI(num), ODF(obj), + fun, flo); + break; + case N_BIGNUM: + ndivide_flonum(num, bi_getd(NRBI(num)), ODF(obj), + fun, flo); + break; + case N_FLONUM: + ndivide_flonum(num, NRFF(num), ODF(obj), fun, flo); + break; + case N_FIXRATIO: + ndivide_flonum(num, + (double)NRFRN(num) / (double)NRFRD(num), + ODF(obj), fun, flo); + break; + case N_BIGRATIO: + ndivide_flonum(num, br_getd(NRBR(num)), ODF(obj), + fun, flo); + break; + } + break; + case LispRatio_t: + switch (NRTYPE(num)) { + case N_FIXNUM: + case N_BIGNUM: + ndivide_xi_xr(num, obj, fun, flo); + break; + case N_FLONUM: + ndivide_flonum(num, NRFF(num), + (double)OFRN(obj) / (double)OFRD(obj), + fun, flo); + break; + case N_FIXRATIO: + case N_BIGRATIO: + ndivide_xr_xr(num, obj, fun, flo); + break; + } + break; + case LispBigratio_t: + switch (NRTYPE(num)) { + case N_FIXNUM: + case N_BIGNUM: + ndivide_xi_xr(num, obj, fun, flo); + break; + case N_FLONUM: + ndivide_flonum(num, NRFF(num), br_getd(OBR(obj)), + fun, flo); + break; + case N_FIXRATIO: + case N_BIGRATIO: + ndivide_xr_xr(num, obj, fun, flo); + break; + } + break; + default: + fatal_object_error(obj, NOT_A_REAL_NUMBER); + break; + } +} + + +/************************************************************************ + * COMPARE + ************************************************************************/ +static int +cmp_real_real(n_real *op1, n_real *op2) +{ + switch (RTYPE(op1)) { + case N_FIXNUM: + switch (RTYPE(op2)) { + case N_FIXNUM: + return (cmp_fi_fi(RFI(op1), RFI(op2))); + case N_BIGNUM: + return (cmp_fi_bi(RFI(op1), RBI(op2))); + case N_FLONUM: + return (cmp_flonum((double)RFI(op1), RFF(op2))); + case N_FIXRATIO: + return (cmp_fi_fr(RFI(op1), RFRN(op2), RFRD(op2))); + case N_BIGRATIO: + return (cmp_fi_br(RFI(op1), RBR(op2))); + } + break; + case N_BIGNUM: + switch (RTYPE(op2)) { + case N_FIXNUM: + return (cmp_bi_fi(RBI(op1), RFI(op2))); + case N_BIGNUM: + return (cmp_bi_bi(RBI(op1), RBI(op2))); + case N_FLONUM: + return (cmp_flonum(bi_getd(RBI(op1)), RFF(op2))); + case N_FIXRATIO: + return (cmp_bi_fr(RBI(op1), RFRN(op2), RFRD(op2))); + case N_BIGRATIO: + return (cmp_bi_br(RBI(op1), RBR(op2))); + } + break; + case N_FLONUM: + switch (RTYPE(op2)) { + case N_FIXNUM: + return (cmp_flonum(RFF(op1), (double)RFI(op2))); + case N_BIGNUM: + return (cmp_flonum(RFF(op1), bi_getd(RBI(op2)))); + case N_FLONUM: + return (cmp_flonum(RFF(op1), RFF(op2))); + case N_FIXRATIO: + return (cmp_flonum(RFF(op1), + (double)RFRN(op2) / (double)RFRD(op2))); + case N_BIGRATIO: + return (cmp_flonum(RFF(op1), br_getd(RBR(op2)))); + } + break; + case N_FIXRATIO: + switch (RTYPE(op2)) { + case N_FIXNUM: + return (cmp_fr_fi(RFRN(op1), RFRD(op1), RFI(op2))); + case N_BIGNUM: + return (cmp_fr_bi(RFRN(op1), RFRD(op1), RBI(op2))); + case N_FLONUM: + return (cmp_flonum((double)RFRN(op1) / (double)RFRD(op1), + RFF(op2))); + case N_FIXRATIO: + return (cmp_fr_fr(RFRN(op1), RFRD(op1), + RFRN(op2), RFRD(op2))); + case N_BIGRATIO: + return (cmp_fr_br(RFRN(op1), RFRD(op1), RBR(op2))); + } + break; + case N_BIGRATIO: + switch (RTYPE(op2)) { + case N_FIXNUM: + return (cmp_br_fi(RBR(op1), RFI(op2))); + case N_BIGNUM: + return (cmp_br_bi(RBR(op1), RBI(op2))); + case N_FLONUM: + return (cmp_flonum(br_getd(RBR(op1)), RFF(op2))); + case N_FIXRATIO: + return (cmp_br_fr(RBR(op1), RFRN(op2), RFRD(op2))); + case N_BIGRATIO: + return (cmp_br_br(RBR(op1), RBR(op2))); + } + } + + return (0); +} + +static int +cmp_real_object(n_real *op1, LispObj *op2) +{ + switch (OBJECT_TYPE(op2)) { + case LispFixnum_t: + switch (RTYPE(op1)) { + case N_FIXNUM: + return (cmp_fi_fi(RFI(op1), OFI(op2))); + case N_BIGNUM: + return (cmp_bi_fi(RBI(op1), OFI(op2))); + case N_FLONUM: + return (cmp_flonum(RFF(op1), (double)OFI(op2))); + case N_FIXRATIO: + return (cmp_fr_fi(RFRD(op1), RFRN(op1), OFI(op2))); + case N_BIGRATIO: + return (cmp_br_fi(RBR(op1), OFI(op2))); + } + break; + case LispInteger_t: + switch (RTYPE(op1)) { + case N_FIXNUM: + return (cmp_fi_fi(RFI(op1), OII(op2))); + case N_BIGNUM: + return (cmp_bi_fi(RBI(op1), OII(op2))); + case N_FLONUM: + return (cmp_flonum(RFF(op1), (double)OII(op2))); + case N_FIXRATIO: + return (cmp_fr_fi(RFRD(op1), RFRN(op1), OII(op2))); + case N_BIGRATIO: + return (cmp_br_fi(RBR(op1), OII(op2))); + } + break; + case LispBignum_t: + switch (RTYPE(op1)) { + case N_FIXNUM: + return (cmp_fi_bi(RFI(op1), OBI(op2))); + case N_BIGNUM: + return (cmp_bi_bi(RBI(op1), OBI(op2))); + case N_FLONUM: + return (cmp_flonum(RFF(op1), bi_getd(OBI(op2)))); + case N_FIXRATIO: + return (cmp_fr_bi(RFRD(op1), RFRN(op1), OBI(op2))); + case N_BIGRATIO: + return (cmp_br_bi(RBR(op1), OBI(op2))); + } + break; + case LispDFloat_t: + switch (RTYPE(op1)) { + case N_FIXNUM: + return (cmp_flonum((double)RFI(op1), ODF(op2))); + case N_BIGNUM: + return (cmp_flonum(bi_getd(RBI(op1)), ODF(op2))); + case N_FLONUM: + return (cmp_flonum(RFF(op1), ODF(op2))); + case N_FIXRATIO: + return (cmp_flonum((double)RFRN(op1) / (double)RFRD(op1), + ODF(op2))); + case N_BIGRATIO: + return (cmp_flonum(br_getd(RBR(op1)), ODF(op2))); + } + break; + case LispRatio_t: + switch (RTYPE(op1)) { + case N_FIXNUM: + return (cmp_fi_fr(RFI(op1), OFRN(op2), OFRD(op2))); + case N_BIGNUM: + return (cmp_bi_fr(RBI(op1), OFRN(op2), OFRD(op2))); + case N_FLONUM: + return (cmp_flonum(RFF(op1), + (double)OFRN(op2) / (double)OFRD(op2))); + case N_FIXRATIO: + return (cmp_fr_fr(RFRN(op1), RFRD(op1), + OFRN(op2), OFRD(op2))); + case N_BIGRATIO: + return (cmp_br_fr(RBR(op1), OFRN(op2), OFRD(op2))); + } + break; + case LispBigratio_t: + switch (RTYPE(op1)) { + case N_FIXNUM: + return (cmp_fi_br(RFI(op1), OBR(op2))); + case N_BIGNUM: + return (cmp_bi_br(RBI(op1), OBR(op2))); + case N_FLONUM: + return (cmp_flonum(RFF(op1), br_getd(OBR(op2)))); + case N_FIXRATIO: + return (cmp_fr_br(RFRN(op1), RFRD(op1), OBR(op2))); + case N_BIGRATIO: + return (cmp_br_br(RBR(op1), OBR(op2))); + } + break; + default: + fatal_object_error(op2, NOT_A_REAL_NUMBER); + break; + } + + return (0); +} + +#if 0 /* not used */ +static int +cmp_number_object(n_number *op1, LispObj *op2) +{ + if (op1->complex) { + if (OBJECT_TYPE(op2) == LispComplex_t) { + if (cmp_real_object(NREAL(op1), OCXR(op2)) == 0) + return (cmp_real_object(NIMAG(op1), OCXI(op2))); + return (1); + } + else if (cmp_real_real(NIMAG(op1), &zero) == 0) + return (cmp_real_object(NREAL(op1), op2)); + else + return (1); + } + else { + switch (OBJECT_TYPE(op2)) { + case LispFixnum_t: + switch (NRTYPE(op1)) { + case N_FIXNUM: + return (cmp_fi_fi(NRFI(op1), OFI(op2))); + case N_BIGNUM: + return (cmp_bi_fi(NRBI(op1), OFI(op2))); + case N_FLONUM: + return (cmp_flonum(NRFF(op1), (double)OFI(op2))); + case N_FIXRATIO: + return (cmp_fr_fi(NRFRD(op1), NRFRN(op1), OFI(op2))); + case N_BIGRATIO: + return (cmp_br_fi(NRBR(op1), OFI(op2))); + } + break; + case LispInteger_t: + switch (NRTYPE(op1)) { + case N_FIXNUM: + return (cmp_fi_fi(NRFI(op1), OII(op2))); + case N_BIGNUM: + return (cmp_bi_fi(NRBI(op1), OII(op2))); + case N_FLONUM: + return (cmp_flonum(NRFF(op1), (double)OII(op2))); + case N_FIXRATIO: + return (cmp_fr_fi(NRFRD(op1), NRFRN(op1), OII(op2))); + case N_BIGRATIO: + return (cmp_br_fi(NRBR(op1), OII(op2))); + } + break; + case LispBignum_t: + switch (NRTYPE(op1)) { + case N_FIXNUM: + return (cmp_fi_bi(NRFI(op1), OBI(op2))); + case N_BIGNUM: + return (cmp_bi_bi(NRBI(op1), OBI(op2))); + case N_FLONUM: + return (cmp_flonum(NRFF(op1), bi_getd(OBI(op2)))); + case N_FIXRATIO: + return (cmp_fr_bi(NRFRD(op1), NRFRN(op1), OBI(op2))); + case N_BIGRATIO: + return (cmp_br_bi(NRBR(op1), OBI(op2))); + } + break; + case LispDFloat_t: + switch (NRTYPE(op1)) { + case N_FIXNUM: + return (cmp_flonum((double)NRFI(op1), ODF(op2))); + case N_BIGNUM: + return (cmp_flonum(bi_getd(NRBI(op1)), ODF(op2))); + case N_FLONUM: + return (cmp_flonum(NRFF(op1), ODF(op2))); + case N_FIXRATIO: + return (cmp_flonum((double)NRFRN(op1) / + (double)NRFRD(op1), + ODF(op2))); + case N_BIGRATIO: + return (cmp_flonum(br_getd(NRBR(op1)), ODF(op2))); + } + break; + case LispRatio_t: + switch (NRTYPE(op1)) { + case N_FIXNUM: + return (cmp_fi_fr(NRFI(op1), OFRN(op2), OFRD(op2))); + case N_BIGNUM: + return (cmp_bi_fr(NRBI(op1), OFRN(op2), OFRD(op2))); + case N_FLONUM: + return (cmp_flonum(NRFF(op1), + (double)OFRN(op2) / (double)OFRD(op2))); + case N_FIXRATIO: + return (cmp_fr_fr(NRFRN(op1), NRFRD(op1), + OFRN(op2), OFRD(op2))); + case N_BIGRATIO: + return (cmp_br_fr(NRBR(op1), OFRN(op2), OFRD(op2))); + } + break; + case LispBigratio_t: + switch (NRTYPE(op1)) { + case N_FIXNUM: + return (cmp_fi_br(NRFI(op1), OBR(op2))); + case N_BIGNUM: + return (cmp_bi_br(NRBI(op1), OBR(op2))); + case N_FLONUM: + return (cmp_flonum(NRFF(op1), br_getd(OBR(op2)))); + case N_FIXRATIO: + return (cmp_fr_br(NRFRN(op1), NRFRD(op1), OBR(op2))); + case N_BIGRATIO: + return (cmp_br_br(NRBR(op1), OBR(op2))); + } + break; + case LispComplex_t: + if (cmp_real_object(&zero, OCXI(op2)) == 0) + return (cmp_real_object(NREAL(op1), OCXR(op2))); + return (1); + default: + fatal_object_error(op2, NOT_A_NUMBER); + break; + } + } + + return (0); +} +#endif + +static int +cmp_object_object(LispObj *op1, LispObj *op2, int real) +{ + if (OBJECT_TYPE(op1) == LispComplex_t) { + if (real) + fatal_object_error(op1, NOT_A_REAL_NUMBER); + if (OBJECT_TYPE(op2) == LispComplex_t) + return (cmp_cx_cx(op1, op2)); + else if (cmp_real_object(&zero, OCXI(op1)) == 0) + return (cmp_object_object(OCXR(op1), op2, real)); + return (1); + } + else if (OBJECT_TYPE(op2) == LispComplex_t) { + if (real) + fatal_object_error(op1, NOT_A_REAL_NUMBER); + if (cmp_real_object(&zero, OCXI(op2)) == 0) + return (cmp_object_object(op1, OCXR(op2), real)); + return (1); + } + else { + switch (OBJECT_TYPE(op1)) { + case LispFixnum_t: + switch (OBJECT_TYPE(op2)) { + case LispFixnum_t: + return (cmp_fi_fi(OFI(op1), OFI(op2))); + case LispInteger_t: + return (cmp_fi_fi(OFI(op1), OII(op2))); + case LispBignum_t: + return (cmp_fi_bi(OFI(op1), OBI(op2))); + case LispDFloat_t: + return (cmp_flonum((double)OFI(op1), ODF(op2))); + case LispRatio_t: + return (cmp_fi_fr(OFI(op1), + OFRN(op2), OFRD(op2))); + case LispBigratio_t: + return (cmp_fi_br(OFI(op1), OBR(op2))); + default: + break; + } + break; + case LispInteger_t: + switch (OBJECT_TYPE(op2)) { + case LispFixnum_t: + return (cmp_fi_fi(OII(op1), OFI(op2))); + case LispInteger_t: + return (cmp_fi_fi(OII(op1), OII(op2))); + case LispBignum_t: + return (cmp_fi_bi(OII(op1), OBI(op2))); + case LispDFloat_t: + return (cmp_flonum((double)OII(op1), ODF(op2))); + case LispRatio_t: + return (cmp_fi_fr(OII(op1), + OFRN(op2), OFRD(op2))); + case LispBigratio_t: + return (cmp_fi_br(OII(op1), OBR(op2))); + default: + break; + } + break; + case LispBignum_t: + switch (OBJECT_TYPE(op2)) { + case LispFixnum_t: + return (cmp_bi_fi(OBI(op1), OFI(op2))); + case LispInteger_t: + return (cmp_bi_fi(OBI(op1), OII(op2))); + case LispBignum_t: + return (cmp_bi_bi(OBI(op1), OBI(op2))); + case LispDFloat_t: + return (cmp_flonum(bi_getd(OBI(op1)), ODF(op2))); + case LispRatio_t: + return (cmp_bi_fr(OBI(op1), + OFRN(op2), OFRD(op2))); + case LispBigratio_t: + return (cmp_bi_br(OBI(op1), OBR(op2))); + default: + break; + } + break; + case LispDFloat_t: + switch (OBJECT_TYPE(op2)) { + case LispFixnum_t: + return (cmp_flonum(ODF(op1), (double)OFI(op2))); + case LispInteger_t: + return (cmp_flonum(ODF(op1), (double)OII(op2))); + case LispBignum_t: + return (cmp_flonum(ODF(op1), bi_getd(OBI(op2)))); + case LispDFloat_t: + return (cmp_flonum(ODF(op1), ODF(op2))); + break; + case LispRatio_t: + return (cmp_flonum(ODF(op1), + (double)OFRN(op2) / + (double)OFRD(op2))); + case LispBigratio_t: + return (cmp_flonum(ODF(op1), br_getd(OBR(op2)))); + default: + break; + } + break; + case LispRatio_t: + switch (OBJECT_TYPE(op2)) { + case LispFixnum_t: + return (cmp_fr_fi(OFRN(op1), OFRD(op1), OFI(op2))); + case LispInteger_t: + return (cmp_fr_fi(OFRN(op1), OFRD(op1), OII(op2))); + case LispBignum_t: + return (cmp_fr_bi(OFRN(op1), OFRD(op1), OBI(op2))); + case LispDFloat_t: + return (cmp_flonum((double)OFRN(op1) / + (double)OFRD(op1), + ODF(op2))); + case LispRatio_t: + return (cmp_fr_fr(OFRN(op1), OFRD(op1), + OFRN(op2), OFRD(op2))); + case LispBigratio_t: + return (cmp_fr_br(OFRN(op1), OFRD(op1), OBR(op2))); + default: + break; + } + break; + case LispBigratio_t: + switch (OBJECT_TYPE(op2)) { + case LispFixnum_t: + return (cmp_br_fi(OBR(op1), OFI(op2))); + case LispInteger_t: + return (cmp_br_fi(OBR(op1), OII(op2))); + case LispBignum_t: + return (cmp_br_bi(OBR(op1), OBI(op2))); + case LispDFloat_t: + return (cmp_flonum(br_getd(OBR(op1)), ODF(op2))); + case LispRatio_t: + return (cmp_br_fr(OBR(op1), OFRN(op2), OFRD(op2))); + case LispBigratio_t: + return (cmp_br_br(OBR(op1), OBR(op2))); + default: + break; + } + break; + default: + fatal_object_error(op1, NOT_A_NUMBER); + break; + } + } + + fatal_object_error(op2, NOT_A_NUMBER); + return (0); +} + + +/************************************************************************ + * FIXNUM + ************************************************************************/ +/* + * check if op1 + op2 will overflow + */ +static INLINE int +fi_fi_add_overflow(long op1, long op2) +{ + long op = op1 + op2; + + return (op1 > 0 ? op2 > op : op2 < op); +} + +/* + * check if op1 - op2 will overflow + */ +static INLINE int +fi_fi_sub_overflow(long op1, long op2) +{ + long op = op1 - op2; + + return (((op1 < 0) ^ (op2 < 0)) && ((op < 0) ^ (op1 < 0))); +} + +/* + * check if op1 * op2 will overflow + */ +static INLINE int +fi_fi_mul_overflow(long op1, long op2) +{ +#ifndef LONG64 + double op = (double)op1 * (double)op2; + + return (op > 2147483647.0 || op < -2147483648.0); +#else + int shift, sign; + long mask; + + if (op1 == 0 || op1 == 1 || op2 == 0 || op2 == 1) + return (0); + + if (op1 == MINSLONG || op2 == MINSLONG) + return (1); + + sign = (op1 < 0) ^ (op2 < 0); + + if (op1 < 0) + op1 = -op1; + if (op2 < 0) + op2 = -op2; + + for (shift = 0, mask = FI_MASK; shift < LONGSBITS; shift++, mask >>= 1) + if (op1 & mask) + break; + ++shift; + for (mask = FI_MASK; shift < LONGSBITS; shift++, mask >>= 1) + if (op2 & mask) + break; + + return (shift < LONGSBITS); +#endif +} + + +/************************************************************************ + * BIGNUM + ************************************************************************/ +static void +rbi_canonicalize(n_real *real) +{ + if (mpi_fiti(RBI(real))) { + long fi = mpi_geti(RBI(real)); + + RTYPE(real) = N_FIXNUM; + mpi_clear(RBI(real)); + XFREE(RBI(real)); + RFI(real) = fi; + } +} + + +/************************************************************************ + * RATIO + ************************************************************************/ +static void +rfr_canonicalize(n_real *real) +{ + long num, numerator, den, denominator, rest; + + num = numerator = RFRN(real); + den = denominator = RFRD(real); + if (denominator == 0) + fatal_error(DIVIDE_BY_ZERO); + + if (num == MINSLONG || den == MINSLONG) { + mpr *bigratio = XALLOC(mpr); + + mpr_init(bigratio); + mpr_seti(bigratio, num, den); + RTYPE(real) = N_BIGRATIO; + RBR(real) = bigratio; + rbr_canonicalize(real); + return; + } + + if (num < 0) + num = -num; + else if (num == 0) { + RFI(real) = 0; + RTYPE(real) = N_FIXNUM; + return; + } + for (;;) { + if ((rest = den % num) == 0) + break; + den = num; + num = rest; + } + if (den != 1) { + denominator /= num; + numerator /= num; + } + if (denominator < 0) { + numerator = -numerator; + denominator = -denominator; + } + if (denominator == 1) { + RTYPE(real) = N_FIXNUM; + RFI(real) = numerator; + } + else { + RFRN(real) = numerator; + RFRD(real) = denominator; + } +} + +static void +rbr_canonicalize(n_real *real) +{ + int fitnum, fitden; + long numerator, denominator; + + mpr_canonicalize(RBR(real)); + fitnum = mpi_fiti(RBRN(real)); + fitden = mpi_fiti(RBRD(real)); + if (fitnum && fitden) { + numerator = mpi_geti(RBRN(real)); + denominator = mpi_geti(RBRD(real)); + mpr_clear(RBR(real)); + XFREE(RBR(real)); + if (numerator == 0) { + RFI(real) = 0; + RTYPE(real) = N_FIXNUM; + } + else if (denominator == 1) { + RTYPE(real) = N_FIXNUM; + RFI(real) = numerator; + } + else { + RTYPE(real) = N_FIXRATIO; + RFRN(real) = numerator; + RFRD(real) = denominator; + } + } + else if (fitden) { + denominator = mpi_geti(RBRD(real)); + if (denominator == 1) { + mpi *bigi = XALLOC(mpi); + + mpi_init(bigi); + mpi_set(bigi, RBRN(real)); + mpr_clear(RBR(real)); + XFREE(RBR(real)); + RTYPE(real) = N_BIGNUM; + RBI(real) = bigi; + } + else if (denominator == 0) + fatal_error(DIVIDE_BY_ZERO); + } +} + + +/************************************************************************ + * COMPLEX + ************************************************************************/ +static void +ncx_canonicalize(n_number *num) +{ + if (NITYPE(num) == N_FIXNUM && NIFI(num) == 0) + num->complex = 0; +} + + +/************************************************************************ + * DIVIDE + ************************************************************************/ +#define NDIVIDE_NOP 0 +#define NDIVIDE_ADD 1 +#define NDIVIDE_SUB 2 +static void +ndivide_fi_fi(n_number *num, long div, int fun, int flo) +{ + long quo, rem; + + if (NRFI(num) == MINSLONG || div == MINSLONG) { + LispObj integer; + mpi *bignum = XALLOC(mpi); + + mpi_init(bignum); + mpi_seti(bignum, NRFI(num)); + NRBI(num) = bignum; + NRTYPE(num) = N_BIGNUM; + integer.type = LispInteger_t; + integer.data.integer = div; + ndivide_xi_xi(num, &integer, fun, flo); + return; + } + else { + quo = NRFI(num) / div; + rem = NRFI(num) % div; + } + + switch (fun) { + case NDIVIDE_CEIL: + if ((rem < 0 && div < 0) || (rem > 0 && div > 0)) { + ++quo; + rem -= div; + } + break; + case NDIVIDE_FLOOR: + if ((rem < 0 && div > 0) || (rem > 0 && div < 0)) { + --quo; + rem += div; + } + break; + case NDIVIDE_ROUND: + if (div > 0) { + if (rem > 0) { + if (rem >= (div + 1) / 2) { + ++quo; + rem -= div; + } + } + else { + if (rem <= (-div - 1) / 2) { + --quo; + rem += div; + } + } + } + else { + if (rem > 0) { + if (rem >= (-div + 1) / 2) { + --quo; + rem += div; + } + } + else { + if (rem <= (div - 1) / 2) { + ++quo; + rem -= div; + } + } + } + break; + } + + NITYPE(num) = N_FIXNUM; + NIFI(num) = rem; + if (flo) { + NRTYPE(num) = N_FLONUM; + NRFF(num) = (double)quo; + } + else + NRFI(num) = quo; +} + +static void +ndivide_xi_xi(n_number *num, LispObj *div, int fun, int flo) +{ + LispType type = OBJECT_TYPE(div); + int state = NDIVIDE_NOP, dsign, rsign; + mpi *quo, *rem; + + quo = XALLOC(mpi); + mpi_init(quo); + if (NRTYPE(num) == N_FIXNUM) + mpi_seti(quo, NRFI(num)); + else + mpi_set(quo, NRBI(num)); + + rem = XALLOC(mpi); + mpi_init(rem); + + switch (type) { + case LispFixnum_t: + mpi_seti(rem, OFI(div)); + break; + case LispInteger_t: + mpi_seti(rem, OII(div)); + break; + default: + mpi_set(rem, OBI(div)); + } + + dsign = mpi_sgn(rem); + + mpi_divqr(quo, rem, quo, rem); + rsign = mpi_sgn(rem); + + switch (fun) { + case NDIVIDE_CEIL: + if ((rsign < 0 && dsign < 0) || (rsign > 0 && dsign > 0)) + state = NDIVIDE_ADD; + break; + case NDIVIDE_FLOOR: + if ((rsign < 0 && dsign > 0) || (rsign > 0 && dsign < 0)) + state = NDIVIDE_SUB; + break; + case NDIVIDE_ROUND: { + mpi test; + + mpi_init(&test); + switch (type) { + case LispFixnum_t: + mpi_seti(&test, OFI(div)); + break; + case LispInteger_t: + mpi_seti(&test, OII(div)); + break; + default: + mpi_set(&test, OBI(div)); + } + if (dsign > 0) { + if (rsign > 0) { + mpi_addi(&test, &test, 1); + mpi_divi(&test, &test, 2); + if (mpi_cmp(rem, &test) >= 0) + state = NDIVIDE_ADD; + } + else { + mpi_neg(&test, &test); + mpi_subi(&test, &test, 1); + mpi_divi(&test, &test, 2); + if (mpi_cmp(rem, &test) <= 0) + state = NDIVIDE_SUB; + } + } + else { + if (rsign > 0) { + mpi_neg(&test, &test); + mpi_addi(&test, &test, 1); + mpi_divi(&test, &test, 2); + if (mpi_cmp(rem, &test) >= 0) + state = NDIVIDE_SUB; + } + else { + mpi_subi(&test, &test, 1); + mpi_divi(&test, &test, 2); + if (mpi_cmp(rem, &test) <= 0) + state = NDIVIDE_ADD; + } + } + mpi_clear(&test); + } break; + } + + if (state == NDIVIDE_ADD) { + mpi_addi(quo, quo, 1); + switch (type) { + case LispFixnum_t: + mpi_subi(rem, rem, OFI(div)); + break; + case LispInteger_t: + mpi_subi(rem, rem, OII(div)); + break; + default: + mpi_sub(rem, rem, OBI(div)); + } + } + else if (state == NDIVIDE_SUB) { + mpi_subi(quo, quo, 1); + switch (type) { + case LispFixnum_t: + mpi_addi(rem, rem, OFI(div)); + break; + case LispInteger_t: + mpi_addi(rem, rem, OII(div)); + break; + default: + mpi_add(rem, rem, OBI(div)); + } + } + + if (mpi_fiti(rem)) { + NITYPE(num) = N_FIXNUM; + NIFI(num) = mpi_geti(rem); + mpi_clear(rem); + XFREE(rem); + } + else { + NITYPE(num) = N_BIGNUM; + NIBI(num) = rem; + } + + clear_real(NREAL(num)); + + if (flo) { + double dval = bi_getd(quo); + + mpi_clear(quo); + XFREE(quo); + NRTYPE(num) = N_FLONUM; + NRFF(num) = dval; + } + else { + NRTYPE(num) = N_BIGNUM; + NRBI(num) = quo; + rbi_canonicalize(NREAL(num)); + } +} + +static void +ndivide_flonum(n_number *number, double num, double div, int fun, int flo) +{ + double quo, rem, modp, tmp; + + modp = modf(num / div, &quo); + rem = num - quo * div; + + switch (fun) { + case NDIVIDE_CEIL: + if ((rem < 0.0 && div < 0.0) || (rem > 0.0 && div > 0.0)) { + quo += 1.0; + rem -= div; + } + break; + case NDIVIDE_FLOOR: + if ((rem < 0.0 && div > 0.0) || (rem > 0.0 && div < 0.0)) { + quo -= 1.0; + rem += div; + } + break; + case NDIVIDE_ROUND: + if (fabs(modp) != 0.5 || modf(quo * 0.5, &tmp) != 0.0) { + if (div > 0.0) { + if (rem > 0.0) { + if (rem >= div * 0.5) { + quo += 1.0; + rem -= div; + } + } + else { + if (rem <= div * -0.5) { + quo -= 1.0; + rem += div; + } + } + } + else { + if (rem > 0.0) { + if (rem >= div * -0.5) { + quo -= 1.0; + rem += div; + } + } + else { + if (rem <= div * 0.5) { + quo += 1.0; + rem -= div; + } + } + } + } + break; + } + if (!finite(quo) || !finite(rem)) + fatal_error(FLOATING_POINT_OVERFLOW); + + NITYPE(number) = N_FLONUM; + NIFF(number) = rem; + + clear_real(NREAL(number)); + + if (flo) { + NRTYPE(number) = N_FLONUM; + NRFF(number) = quo; + } + else { + if ((long)quo == quo) { + NRTYPE(number) = N_FIXNUM; + NRFI(number) = (long)quo; + } + else { + mpi *bigi = XALLOC(mpi); + + mpi_init(bigi); + mpi_setd(bigi, quo); + NRBI(number) = bigi; + NRTYPE(number) = N_BIGNUM; + } + } +} + +static void +ndivide_xi_xr(n_number *num, LispObj *div, int fun, int flo) +{ + int state = NDIVIDE_NOP, dsign, rsign; + mpi *quo; + mpr *rem; + + quo = XALLOC(mpi); + mpi_init(quo); + if (NRTYPE(num) == N_FIXNUM) + mpi_seti(quo, NRFI(num)); + else + mpi_set(quo, NRBI(num)); + + rem = XALLOC(mpr); + mpr_init(rem); + + if (XOBJECT_TYPE(div) == LispRatio_t) + mpr_seti(rem, OFRN(div), OFRD(div)); + else + mpr_set(rem, OBR(div)); + dsign = mpi_sgn(mpr_num(rem)); + mpi_mul(quo, quo, mpr_den(rem)); + + mpi_divqr(quo, mpr_num(rem), quo, mpr_num(rem)); + mpr_canonicalize(rem); + + rsign = mpi_sgn(mpr_num(rem)); + if (mpr_fiti(rem)) { + if (mpi_geti(mpr_den(rem)) == 1) { + NITYPE(num) = N_FIXNUM; + NIFI(num) = mpi_geti(mpr_num(rem)); + } + else { + NITYPE(num) = N_FIXRATIO; + NIFRN(num) = mpi_geti(mpr_num(rem)); + NIFRD(num) = mpi_geti(mpr_den(rem)); + } + mpr_clear(rem); + XFREE(rem); + } + else { + if (mpi_fiti(mpr_den(rem)) && mpi_geti(mpr_den(rem)) == 1) { + NITYPE(num) = N_BIGNUM; + NIBI(num) = mpr_num(rem); + mpi_clear(mpr_den(rem)); + XFREE(rem); + } + else { + NITYPE(num) = N_BIGRATIO; + NIBR(num) = rem; + } + } + + switch (fun) { + case NDIVIDE_CEIL: + if ((rsign < 0 && dsign < 0) || (rsign > 0 && dsign > 0)) + state = NDIVIDE_ADD; + break; + case NDIVIDE_FLOOR: + if ((rsign < 0 && dsign > 0) || (rsign > 0 && dsign < 0)) + state = NDIVIDE_SUB; + break; + case NDIVIDE_ROUND: { + n_real cmp; + + set_real_object(&cmp, div); + div_real_real(&cmp, &two); + if (dsign > 0) { + if (rsign > 0) { + if (cmp_real_real(NIMAG(num), &cmp) >= 0) + state = NDIVIDE_ADD; + } + else { + neg_real(&cmp); + if (cmp_real_real(NIMAG(num), &cmp) <= 0) + state = NDIVIDE_SUB; + } + } + else { + if (rsign > 0) { + neg_real(&cmp); + if (cmp_real_real(NIMAG(num), &cmp) >= 0) + state = NDIVIDE_SUB; + } + else { + if (cmp_real_real(NIMAG(num), &cmp) <= 0) + state = NDIVIDE_ADD; + } + } + clear_real(&cmp); + } break; + } + + if (state == NDIVIDE_ADD) { + mpi_addi(quo, quo, 1); + sub_real_object(NIMAG(num), div); + } + else if (state == NDIVIDE_SUB) { + mpi_subi(quo, quo, 1); + add_real_object(NIMAG(num), div); + } + + clear_real(NREAL(num)); + + if (flo) { + double dval = bi_getd(quo); + + mpi_clear(quo); + XFREE(quo); + NRTYPE(num) = N_FLONUM; + NRFF(num) = dval; + } + else { + NRBI(num) = quo; + NRTYPE(num) = N_BIGNUM; + rbi_canonicalize(NREAL(num)); + } +} + +static void +ndivide_xr_xi(n_number *num, LispObj *div, int fun, int flo) +{ + LispType type = OBJECT_TYPE(div); + int state = NDIVIDE_NOP, dsign, rsign; + mpi *quo; + mpr *rem; + + quo = XALLOC(mpi); + mpi_init(quo); + switch (type) { + case LispFixnum_t: + dsign = OFI(div) < 0 ? -1 : OFI(div) > 0 ? 1 : 0; + mpi_seti(quo, OFI(div)); + break; + case LispInteger_t: + dsign = OII(div) < 0 ? -1 : OII(div) > 0 ? 1 : 0; + mpi_seti(quo, OII(div)); + break; + default: + dsign = mpi_sgn(OBI(div)); + mpi_set(quo, OBI(div)); + break; + } + + rem = XALLOC(mpr); + mpr_init(rem); + if (NRTYPE(num) == N_FIXRATIO) { + mpr_seti(rem, NRFRN(num), NRFRD(num)); + mpi_muli(quo, quo, NRFRD(num)); + } + else { + mpr_set(rem, NRBR(num)); + mpi_mul(quo, quo, NRBRD(num)); + } + mpi_divqr(quo, mpr_num(rem), mpr_num(rem), quo); + mpr_canonicalize(rem); + + rsign = mpi_sgn(mpr_num(rem)); + if (mpr_fiti(rem)) { + NITYPE(num) = N_FIXRATIO; + NIFRN(num) = mpi_geti(mpr_num(rem)); + NIFRD(num) = mpi_geti(mpr_den(rem)); + mpr_clear(rem); + XFREE(rem); + } + else { + NITYPE(num) = N_BIGRATIO; + NIBR(num) = rem; + } + + switch (fun) { + case NDIVIDE_CEIL: + if ((rsign < 0 && dsign < 0) || (rsign > 0 && dsign > 0)) + state = NDIVIDE_ADD; + break; + case NDIVIDE_FLOOR: + if ((rsign < 0 && dsign > 0) || (rsign > 0 && dsign < 0)) + state = NDIVIDE_SUB; + break; + case NDIVIDE_ROUND: { + n_real cmp; + + set_real_object(&cmp, div); + div_real_real(&cmp, &two); + if (dsign > 0) { + if (rsign > 0) { + if (cmp_real_real(NIMAG(num), &cmp) >= 0) + state = NDIVIDE_ADD; + } + else { + neg_real(&cmp); + if (cmp_real_real(NIMAG(num), &cmp) <= 0) + state = NDIVIDE_SUB; + } + } + else { + if (rsign > 0) { + neg_real(&cmp); + if (cmp_real_real(NIMAG(num), &cmp) >= 0) + state = NDIVIDE_SUB; + } + else { + if (cmp_real_real(NIMAG(num), &cmp) <= 0) + state = NDIVIDE_ADD; + } + } + clear_real(&cmp); + } break; + } + + if (state == NDIVIDE_ADD) { + mpi_addi(quo, quo, 1); + sub_real_object(NIMAG(num), div); + } + else if (state == NDIVIDE_SUB) { + mpi_subi(quo, quo, 1); + add_real_object(NIMAG(num), div); + } + + clear_real(NREAL(num)); + + if (flo) { + double dval = bi_getd(quo); + + mpi_clear(quo); + XFREE(quo); + NRTYPE(num) = N_FLONUM; + NRFF(num) = dval; + } + else { + NRBI(num) = quo; + NRTYPE(num) = N_BIGNUM; + rbi_canonicalize(NREAL(num)); + } +} + +static void +ndivide_xr_xr(n_number *num, LispObj *div, int fun, int flo) +{ + int state = NDIVIDE_NOP, dsign, rsign, modp; + mpr *bigr; + mpi *bigi; + + bigr = XALLOC(mpr); + mpr_init(bigr); + if (NRTYPE(num) == N_FIXRATIO) + mpr_seti(bigr, NRFRN(num), NRFRD(num)); + else + mpr_set(bigr, NRBR(num)); + + NITYPE(num) = N_BIGRATIO; + NIBR(num) = bigr; + + if (OBJECT_TYPE(div) == LispRatio_t) { + dsign = OFRN(div) < 0 ? -1 : OFRN(div) > 0 ? 1 : 0; + mpi_muli(mpr_num(bigr), mpr_num(bigr), OFRD(div)); + mpi_muli(mpr_den(bigr), mpr_den(bigr), OFRN(div)); + } + else { + dsign = mpi_sgn(OBRN(div)); + mpr_div(bigr, bigr, OBR(div)); + } + modp = mpi_fiti(mpr_den(bigr)) && mpi_geti(mpr_den(bigr)) == 2; + + bigi = XALLOC(mpi); + mpi_init(bigi); + mpi_divqr(bigi, mpr_num(bigr), mpr_num(bigr), mpr_den(bigr)); + + if (OBJECT_TYPE(div) == LispRatio_t) + mpi_seti(mpr_den(bigr), OFRD(div)); + else + mpi_set(mpr_den(bigr), OBRD(div)); + if (NRTYPE(num) == N_FIXRATIO) + mpi_muli(mpr_den(bigr), mpr_den(bigr), NRFRD(num)); + else + mpi_mul(mpr_den(bigr), mpr_den(bigr), NRBRD(num)); + + clear_real(NREAL(num)); + NRTYPE(num) = N_BIGNUM; + NRBI(num) = bigi; + + rbr_canonicalize(NIMAG(num)); + rsign = cmp_real_real(NIMAG(num), &zero); + + switch (fun) { + case NDIVIDE_CEIL: + if ((rsign < 0 && dsign < 0) || (rsign > 0 && dsign > 0)) + state = NDIVIDE_ADD; + break; + case NDIVIDE_FLOOR: + if ((rsign < 0 && dsign > 0) || (rsign > 0 && dsign < 0)) + state = NDIVIDE_SUB; + break; + case NDIVIDE_ROUND: + if (!modp || (bigi->digs[0] & 1) == 1) { + n_real cmp; + + set_real_object(&cmp, div); + div_real_real(&cmp, &two); + if (dsign > 0) { + if (rsign > 0) { + if (cmp_real_real(NIMAG(num), &cmp) >= 0) + state = NDIVIDE_ADD; + } + else { + neg_real(&cmp); + if (cmp_real_real(NIMAG(num), &cmp) <= 0) + state = NDIVIDE_SUB; + } + } + else { + if (rsign > 0) { + neg_real(&cmp); + if (cmp_real_real(NIMAG(num), &cmp) >= 0) + state = NDIVIDE_SUB; + } + else { + if (cmp_real_real(NIMAG(num), &cmp) <= 0) + state = NDIVIDE_ADD; + } + } + clear_real(&cmp); + } + break; + } + + if (state == NDIVIDE_ADD) { + add_real_real(NREAL(num), &one); + sub_real_object(NIMAG(num), div); + } + else if (state == NDIVIDE_SUB) { + sub_real_real(NREAL(num), &one); + add_real_object(NIMAG(num), div); + } + + if (NRTYPE(num) == N_BIGNUM) { + if (flo) { + double dval = bi_getd(bigi); + + mpi_clear(bigi); + XFREE(bigi); + NRTYPE(num) = N_FLONUM; + NRFF(num) = dval; + } + else + rbi_canonicalize(NREAL(num)); + } + else if (flo) { + NRTYPE(num) = N_FLONUM; + NRFF(num) = (double)NRFI(num); + } +} + + +/************************************************************************ + * REAL COMPLEX + ************************************************************************/ +static void +nadd_re_cx(n_number *num, LispObj *comp) +{ +/* + Ra+Rb Ib + */ + /* Ra+Rb */ + add_real_object(NREAL(num), OCXR(comp)); + + /* Ib */ + set_real_object(NIMAG(num), OCXI(comp)); + + num->complex = 1; + + ncx_canonicalize(num); +} + +static void +nsub_re_cx(n_number *num, LispObj *comp) +{ +/* + Ra-Rb -Ib + */ + /* Ra-Rb */ + sub_real_object(NREAL(num), OCXR(comp)); + + /* -Ib */ + NITYPE(num) = N_FIXNUM; + NIFI(num) = -1; + mul_real_object(NIMAG(num), OCXI(comp)); + + num->complex = 1; + + ncx_canonicalize(num); +} + +static void +nmul_re_cx(n_number *num, LispObj *comp) +{ +/* + Ra*Rb Ra*Ib + */ + /* copy before change */ + set_real_real(NIMAG(num), NREAL(num)); + + /* Ra*Rb */ + mul_real_object(NREAL(num), OCXR(comp)); + + /* Ra*Ib */ + mul_real_object(NIMAG(num), OCXI(comp)); + + num->complex = 1; + + ncx_canonicalize(num); +} + +static void +ndiv_re_cx(n_number *num, LispObj *comp) +{ +/* + Ra*Rb -Ib*Ra + ----------- ----------- + Rb*Rb+Ib*Ib Rb*Rb+Ib*Ib + */ + n_real div, temp; + + /* Rb*Rb */ + set_real_object(&div, OCXR(comp)); + mul_real_object(&div, OCXR(comp)); + + /* Ib*Ib */ + set_real_object(&temp, OCXI(comp)); + mul_real_object(&temp, OCXI(comp)); + + /* Rb*Rb+Ib*Ib */ + add_real_real(&div, &temp); + clear_real(&temp); + + /* -Ib*Ra */ + NITYPE(num) = N_FIXNUM; + NIFI(num) = -1; + mul_real_object(NIMAG(num), OCXI(comp)); + mul_real_real(NIMAG(num), NREAL(num)); + + /* Ra*Rb */ + mul_real_object(NREAL(num), OCXR(comp)); + + div_real_real(NREAL(num), &div); + div_real_real(NIMAG(num), &div); + clear_real(&div); + + num->complex = 1; + + ncx_canonicalize(num); +} + + +/************************************************************************ + * COMPLEX REAL + ************************************************************************/ +static void +nadd_cx_re(n_number *num, LispObj *re) +{ +/* + Ra+Rb Ia + */ + add_real_object(NREAL(num), re); + + ncx_canonicalize(num); +} + +static void +nsub_cx_re(n_number *num, LispObj *re) +{ +/* + Ra-Rb Ia + */ + sub_real_object(NREAL(num), re); + + ncx_canonicalize(num); +} + +static void +nmul_cx_re(n_number *num, LispObj *re) +{ +/* + Ra*Rb Ia*Rb + */ + mul_real_object(NREAL(num), re); + mul_real_object(NIMAG(num), re); + + ncx_canonicalize(num); +} + +static void +ndiv_cx_re(n_number *num, LispObj *re) +{ +/* + Ra/Rb Ia/Rb + */ + div_real_object(NREAL(num), re); + div_real_object(NIMAG(num), re); + + ncx_canonicalize(num); +} + + +/************************************************************************ + * COMPLEX COMPLEX + ************************************************************************/ +static void +nadd_cx_cx(n_number *num, LispObj *comp) +{ +/* + Ra+Rb Ia+Ib + */ + add_real_object(NREAL(num), OCXR(comp)); + add_real_object(NIMAG(num), OCXI(comp)); + + ncx_canonicalize(num); +} + +static void +nsub_cx_cx(n_number *num, LispObj *comp) +{ +/* + Ra-Rb Ia-Ib + */ + sub_real_object(NREAL(num), OCXR(comp)); + sub_real_object(NIMAG(num), OCXI(comp)); + + ncx_canonicalize(num); +} + +static void +nmul_cx_cx(n_number *num, LispObj *comp) +{ +/* + Ra*Rb-Ia*Ib Ra*Ib+Ia*Rb + */ + n_real IaIb, RaIb; + + set_real_real(&IaIb, NIMAG(num)); + mul_real_object(&IaIb, OCXI(comp)); + + set_real_real(&RaIb, NREAL(num)); + mul_real_object(&RaIb, OCXI(comp)); + + /* Ra*Rb-Ia*Ib */ + mul_real_object(NREAL(num), OCXR(comp)); + sub_real_real(NREAL(num), &IaIb); + clear_real(&IaIb); + + /* Ra*Ib+Ia*Rb */ + mul_real_object(NIMAG(num), OCXR(comp)); + add_real_real(NIMAG(num), &RaIb); + clear_real(&RaIb); + + ncx_canonicalize(num); +} + +static void +ndiv_cx_cx(n_number *num, LispObj *comp) +{ +/* + Ra*Rb+Ia*Ib Ia*Rb-Ib*Ra + ----------- ----------- + Rb*Rb+Ib*Ib Rb*Rb+Ib*Ib + */ + n_real temp1, temp2; + + /* IaIb */ + set_real_real(&temp1, NIMAG(num)); + mul_real_object(&temp1, OCXI(comp)); + + /* IbRa */ + set_real_real(&temp2, NREAL(num)); + mul_real_object(&temp2, OCXI(comp)); + + /* Ra*Rb+Ia*Ib */ + mul_real_object(NREAL(num), OCXR(comp)); + add_real_real(NREAL(num), &temp1); + clear_real(&temp1); + + /* Ia*Rb-Ib*Ra */ + mul_real_object(NIMAG(num), OCXR(comp)); + sub_real_real(NIMAG(num), &temp2); + clear_real(&temp2); + + + /* Rb*Rb */ + set_real_object(&temp1, OCXR(comp)); + mul_real_object(&temp1, OCXR(comp)); + + /* Ib*Ib */ + set_real_object(&temp2, OCXI(comp)); + mul_real_object(&temp2, OCXI(comp)); + + /* Rb*Rb+Ib*Ib */ + add_real_real(&temp1, &temp2); + clear_real(&temp2); + + div_real_real(NREAL(num), &temp1); + div_real_real(NIMAG(num), &temp1); + clear_real(&temp1); + + ncx_canonicalize(num); +} + +static int +cmp_cx_cx(LispObj *op1, LispObj *op2) +{ + int cmp; + + cmp = cmp_object_object(OCXR(op1), OCXR(op2), 1); + if (cmp == 0) + cmp = cmp_object_object(OCXI(op1), OCXI(op2), 1); + + return (cmp); +} + + +/************************************************************************ + * FLONUM FLONUM + ************************************************************************/ +static void +radd_flonum(n_real *real, double op1, double op2) +{ + double value = op1 + op2; + + if (!finite(value)) + fatal_error(FLOATING_POINT_OVERFLOW); + switch (RTYPE(real)) { + case N_FIXNUM: + case N_FIXRATIO: + RTYPE(real) = N_FLONUM; + break; + case N_BIGNUM: + RCLEAR_BI(real); + RTYPE(real) = N_FLONUM; + break; + case N_BIGRATIO: + RCLEAR_BR(real); + RTYPE(real) = N_FLONUM; + break; + } + RFF(real) = value; +} + +static void +rsub_flonum(n_real *real, double op1, double op2) +{ + double value = op1 - op2; + + if (!finite(value)) + fatal_error(FLOATING_POINT_OVERFLOW); + switch (RTYPE(real)) { + case N_FIXNUM: + case N_FIXRATIO: + RTYPE(real) = N_FLONUM; + break; + case N_BIGNUM: + RCLEAR_BI(real); + RTYPE(real) = N_FLONUM; + break; + case N_BIGRATIO: + RCLEAR_BR(real); + RTYPE(real) = N_FLONUM; + break; + } + RFF(real) = value; +} + +static void +rmul_flonum(n_real *real, double op1, double op2) +{ + double value = op1 * op2; + + if (!finite(value)) + fatal_error(FLOATING_POINT_OVERFLOW); + switch (RTYPE(real)) { + case N_FIXNUM: + case N_FIXRATIO: + RTYPE(real) = N_FLONUM; + break; + case N_BIGNUM: + RCLEAR_BI(real); + RTYPE(real) = N_FLONUM; + break; + case N_BIGRATIO: + RCLEAR_BR(real); + RTYPE(real) = N_FLONUM; + break; + } + RFF(real) = value; +} + +static void +rdiv_flonum(n_real *real, double op1, double op2) +{ + double value; + + if (op2 == 0.0) + fatal_error(DIVIDE_BY_ZERO); + value = op1 / op2; + if (!finite(value)) + fatal_error(FLOATING_POINT_OVERFLOW); + switch (RTYPE(real)) { + case N_FIXNUM: + case N_FIXRATIO: + RTYPE(real) = N_FLONUM; + break; + case N_BIGNUM: + RCLEAR_BI(real); + RTYPE(real) = N_FLONUM; + break; + case N_BIGRATIO: + RCLEAR_BR(real); + RTYPE(real) = N_FLONUM; + break; + } + RFF(real) = value; +} + +static int +cmp_flonum(double op1, double op2) +{ + double value = op1 - op2; + + if (!finite(value)) + fatal_error(FLOATING_POINT_OVERFLOW); + + return (value > 0.0 ? 1 : value < 0.0 ? -1 : 0); +} + + +/************************************************************************ + * FIXNUM FIXNUM + ************************************************************************/ +static void +rop_fi_fi_bi(n_real *real, long fi, int op) +{ + mpi *bigi = XALLOC(mpi); + + mpi_init(bigi); + mpi_seti(bigi, RFI(real)); + if (op == NOP_ADD) + mpi_addi(bigi, bigi, fi); + else if (op == NOP_SUB) + mpi_subi(bigi, bigi, fi); + else + mpi_muli(bigi, bigi, fi); + RBI(real) = bigi; + RTYPE(real) = N_BIGNUM; +} + +static INLINE void +radd_fi_fi(n_real *real, long fi) +{ + if (!fi_fi_add_overflow(RFI(real), fi)) + RFI(real) += fi; + else + rop_fi_fi_bi(real, fi, NOP_ADD); +} + +static INLINE void +rsub_fi_fi(n_real *real, long fi) +{ + if (!fi_fi_sub_overflow(RFI(real), fi)) + RFI(real) -= fi; + else + rop_fi_fi_bi(real, fi, NOP_SUB); +} + +static INLINE void +rmul_fi_fi(n_real *real, long fi) +{ + if (!fi_fi_mul_overflow(RFI(real), fi)) + RFI(real) *= fi; + else + rop_fi_fi_bi(real, fi, NOP_MUL); +} + +static INLINE void +rdiv_fi_fi(n_real *real, long fi) +{ + RTYPE(real) = N_FIXRATIO; + RFRN(real) = RFI(real); + RFRD(real) = fi; + rfr_canonicalize(real); +} + +static INLINE int +cmp_fi_fi(long op1, long op2) +{ + if (op1 > op2) + return (1); + else if (op1 < op2) + return (-1); + + return (0); +} + + +/************************************************************************ + * FIXNUM BIGNUM + ************************************************************************/ +static void +rop_fi_bi_xi(n_real *real, mpi *bi, int nop) +{ + mpi *bigi = XALLOC(mpi); + + mpi_init(bigi); + mpi_seti(bigi, RFI(real)); + if (nop == NOP_ADD) + mpi_add(bigi, bigi, bi); + else if (nop == NOP_SUB) + mpi_sub(bigi, bigi, bi); + else + mpi_mul(bigi, bigi, bi); + + if (mpi_fiti(bigi)) { + RFI(real) = mpi_geti(bigi); + mpi_clear(bigi); + XFREE(bigi); + } + else { + RBI(real) = bigi; + RTYPE(real) = N_BIGNUM; + } +} + +static INLINE void +radd_fi_bi(n_real *real, mpi *bi) +{ + rop_fi_bi_xi(real, bi, NOP_ADD); +} + +static INLINE void +rsub_fi_bi(n_real *real, mpi *bi) +{ + rop_fi_bi_xi(real, bi, NOP_SUB); +} + +static INLINE void +rmul_fi_bi(n_real *real, mpi *bi) +{ + rop_fi_bi_xi(real, bi, NOP_MUL); +} + +static void +rdiv_fi_bi(n_real *real, mpi *bi) +{ + mpr *bigr; + + if (mpi_cmpi(bi, 0) == 0) + fatal_error(DIVIDE_BY_ZERO); + + bigr = XALLOC(mpr); + mpr_init(bigr); + mpi_seti(mpr_num(bigr), RFI(real)); + mpi_set(mpr_den(bigr), bi); + RBR(real) = bigr; + RTYPE(real) = N_BIGRATIO; + rbr_canonicalize(real); +} + +static INLINE int +cmp_fi_bi(long fixnum, mpi *bignum) +{ + return (-mpi_cmpi(bignum, fixnum)); +} + + +/************************************************************************ + * FIXNUM FIXRATIO + ************************************************************************/ +static void +rop_fi_fr_as_xr(n_real *real, long num, long den, int nop) +{ + int fit; + long value = 0, op = RFI(real); + + fit = !fi_fi_mul_overflow(op, den); + if (fit) { + value = op * den; + if (nop == NOP_ADD) + fit = !fi_fi_add_overflow(value, num); + else + fit = !fi_fi_sub_overflow(value, num); + } + if (fit) { + if (nop == NOP_ADD) + RFRN(real) = value + num; + else + RFRN(real) = value - num; + RFRD(real) = den; + RTYPE(real) = N_FIXRATIO; + rfr_canonicalize(real); + } + else { + mpi iop; + mpr *bigr = XALLOC(mpr); + + mpi_init(&iop); + mpi_seti(&iop, op); + mpi_muli(&iop, &iop, den); + + mpr_init(bigr); + mpr_seti(bigr, num, den); + if (nop == NOP_ADD) + mpi_add(mpr_num(bigr), &iop, mpr_num(bigr)); + else + mpi_sub(mpr_num(bigr), &iop, mpr_num(bigr)); + mpi_clear(&iop); + RBR(real) = bigr; + RTYPE(real) = N_BIGRATIO; + rbr_canonicalize(real); + } +} + +static void +rop_fi_fr_md_xr(n_real *real, long num, long den, int nop) +{ + int fit; + long op = RFI(real); + + if (nop == NOP_MUL) + fit = !fi_fi_mul_overflow(op, num); + else + fit = !fi_fi_mul_overflow(op, den); + if (fit) { + if (nop == NOP_MUL) { + RFRN(real) = op * num; + RFRD(real) = den; + } + else { + RFRN(real) = op * den; + RFRD(real) = num; + } + RTYPE(real) = N_FIXRATIO; + rfr_canonicalize(real); + } + else { + mpi iop; + mpr *bigr = XALLOC(mpr); + + mpi_init(&iop); + mpi_seti(&iop, op); + + mpr_init(bigr); + if (nop == NOP_MUL) + mpr_seti(bigr, num, den); + else + mpr_seti(bigr, den, num); + mpi_mul(mpr_num(bigr), mpr_num(bigr), &iop); + mpi_clear(&iop); + RBR(real) = bigr; + RTYPE(real) = N_BIGRATIO; + rbr_canonicalize(real); + } +} + +static INLINE void +radd_fi_fr(n_real *real, long num, long den) +{ + rop_fi_fr_as_xr(real, num, den, NOP_ADD); +} + +static INLINE void +rsub_fi_fr(n_real *real, long num, long den) +{ + rop_fi_fr_as_xr(real, num, den, NOP_SUB); +} + +static INLINE void +rmul_fi_fr(n_real *real, long num, long den) +{ + rop_fi_fr_md_xr(real, num, den, NOP_MUL); +} + +static INLINE void +rdiv_fi_fr(n_real *real, long num, long den) +{ + rop_fi_fr_md_xr(real, num, den, NOP_DIV); +} + +static INLINE int +cmp_fi_fr(long fi, long num, long den) +{ + return (cmp_flonum((double)fi, (double)num / (double)den)); +} + + +/************************************************************************ + * FIXNUM BIGRATIO + ************************************************************************/ +static void +rop_fi_br_as_xr(n_real *real, mpr *ratio, int nop) +{ + mpi iop; + mpr *bigr = XALLOC(mpr); + + mpi_init(&iop); + mpi_seti(&iop, RFI(real)); + + mpr_init(bigr); + mpr_set(bigr, ratio); + + mpi_mul(&iop, &iop, mpr_den(ratio)); + if (nop == NOP_ADD) + mpi_add(mpr_num(bigr), &iop, mpr_num(bigr)); + else + mpi_sub(mpr_num(bigr), &iop, mpr_num(bigr)); + + mpi_clear(&iop); + RBR(real) = bigr; + RTYPE(real) = N_BIGRATIO; + rbr_canonicalize(real); +} + +static void +rop_fi_br_md_xr(n_real *real, mpr *ratio, int nop) +{ + mpi iop; + mpr *bigr = XALLOC(mpr); + + mpi_init(&iop); + mpi_seti(&iop, RFI(real)); + + mpr_init(bigr); + if (nop == NOP_MUL) + mpr_set(bigr, ratio); + else + mpr_inv(bigr, ratio); + + mpi_mul(mpr_num(bigr), &iop, mpr_num(bigr)); + + mpi_clear(&iop); + RBR(real) = bigr; + RTYPE(real) = N_BIGRATIO; + rbr_canonicalize(real); +} + +static INLINE void +radd_fi_br(n_real *real, mpr *ratio) +{ + rop_fi_br_as_xr(real, ratio, NOP_ADD); +} + +static INLINE void +rsub_fi_br(n_real *real, mpr *ratio) +{ + rop_fi_br_as_xr(real, ratio, NOP_SUB); +} + +static INLINE void +rmul_fi_br(n_real *real, mpr *ratio) +{ + rop_fi_br_md_xr(real, ratio, NOP_MUL); +} + +static INLINE void +rdiv_fi_br(n_real *real, mpr *ratio) +{ + rop_fi_br_md_xr(real, ratio, NOP_DIV); +} + +static INLINE int +cmp_fi_br(long op1, mpr *op2) +{ + return (-mpr_cmpi(op2, op1)); +} + + +/************************************************************************ + * BIGNUM FIXNUM + ************************************************************************/ +static INLINE void +radd_bi_fi(n_real *real, long fi) +{ + mpi_addi(RBI(real), RBI(real), fi); + rbi_canonicalize(real); +} + +static INLINE void +rsub_bi_fi(n_real *real, long fi) +{ + mpi_subi(RBI(real), RBI(real), fi); + rbi_canonicalize(real); +} + +static INLINE void +rmul_bi_fi(n_real *real, long fi) +{ + mpi_muli(RBI(real), RBI(real), fi); + rbi_canonicalize(real); +} + +static void +rdiv_bi_fi(n_real *real, long fi) +{ + mpr *bigr; + + if (RFI(real) == 0) + fatal_error(DIVIDE_BY_ZERO); + + bigr = XALLOC(mpr); + mpr_init(bigr); + mpi_set(mpr_num(bigr), RBI(real)); + mpi_seti(mpr_den(bigr), fi); + RCLEAR_BI(real); + RBR(real) = bigr; + RTYPE(real) = N_BIGRATIO; + rbr_canonicalize(real); +} + +static INLINE int +cmp_bi_fi(mpi *bignum, long fi) +{ + return (mpi_cmpi(bignum, fi)); +} + + +/************************************************************************ + * BIGNUM BIGNUM + ************************************************************************/ +static INLINE void +radd_bi_bi(n_real *real, mpi *bignum) +{ + mpi_add(RBI(real), RBI(real), bignum); + rbi_canonicalize(real); +} + +static INLINE void +rsub_bi_bi(n_real *real, mpi *bignum) +{ + mpi_sub(RBI(real), RBI(real), bignum); + rbi_canonicalize(real); +} + +static INLINE void +rmul_bi_bi(n_real *real, mpi *bignum) +{ + mpi_mul(RBI(real), RBI(real), bignum); + rbi_canonicalize(real); +} + +static void +rdiv_bi_bi(n_real *real, mpi *bignum) +{ + mpr *bigr; + + if (mpi_cmpi(bignum, 0) == 0) + fatal_error(DIVIDE_BY_ZERO); + + bigr = XALLOC(mpr); + mpr_init(bigr); + mpi_set(mpr_num(bigr), RBI(real)); + mpi_set(mpr_den(bigr), bignum); + RCLEAR_BI(real); + RBR(real) = bigr; + RTYPE(real) = N_BIGRATIO; + rbr_canonicalize(real); +} + +static INLINE int +cmp_bi_bi(mpi *op1, mpi *op2) +{ + return (mpi_cmp(op1, op2)); +} + + +/************************************************************************ + * BIGNUM FIXRATIO + ************************************************************************/ +static void +rop_bi_fr_as_xr(n_real *real, long num, long den, int nop) +{ + mpi iop; + mpr *bigr = XALLOC(mpr); + + mpi_init(&iop); + mpi_set(&iop, RBI(real)); + mpi_muli(&iop, &iop, den); + + mpr_init(bigr); + mpr_seti(bigr, num, den); + + if (nop == NOP_ADD) + mpi_add(mpr_num(bigr), &iop, mpr_num(bigr)); + else + mpi_sub(mpr_num(bigr), &iop, mpr_num(bigr)); + mpi_clear(&iop); + + RCLEAR_BI(real); + RBR(real) = bigr; + RTYPE(real) = N_BIGRATIO; + rbr_canonicalize(real); +} + +static INLINE void +rop_bi_fr_md_xr(n_real *real, long num, long den, int nop) +{ + mpr *bigr = XALLOC(mpr); + + mpr_init(bigr); + + mpr_seti(bigr, num, den); + + if (nop == NOP_MUL) + mpi_mul(mpr_num(bigr), RBI(real), mpr_num(bigr)); + else { + mpi_mul(mpr_den(bigr), RBI(real), mpr_den(bigr)); + mpr_inv(bigr, bigr); + } + + RCLEAR_BI(real); + RBR(real) = bigr; + RTYPE(real) = N_BIGRATIO; + rbr_canonicalize(real); +} + +static INLINE void +radd_bi_fr(n_real *real, long num, long den) +{ + rop_bi_fr_as_xr(real, num, den, NOP_ADD); +} + +static INLINE void +rsub_bi_fr(n_real *real, long num, long den) +{ + rop_bi_fr_as_xr(real, num, den, NOP_SUB); +} + +static INLINE void +rmul_bi_fr(n_real *real, long num, long den) +{ + rop_bi_fr_md_xr(real, num, den, NOP_MUL); +} + +static INLINE void +rdiv_bi_fr(n_real *real, long num, long den) +{ + rop_bi_fr_md_xr(real, num, den, NOP_DIV); +} + +static int +cmp_bi_fr(mpi *bignum, long num, long den) +{ + int cmp; + mpr cmp1, cmp2; + + mpr_init(&cmp1); + mpi_set(mpr_num(&cmp1), bignum); + mpi_seti(mpr_den(&cmp1), 1); + + mpr_init(&cmp2); + mpr_seti(&cmp2, num, den); + + cmp = mpr_cmp(&cmp1, &cmp2); + mpr_clear(&cmp1); + mpr_clear(&cmp2); + + return (cmp); +} + + +/************************************************************************ + * BIGNUM BIGRATIO + ************************************************************************/ +static void +rop_bi_br_as_xr(n_real *real, mpr *bigratio, int nop) +{ + mpi iop; + mpr *bigr = XALLOC(mpr); + + mpi_init(&iop); + mpi_set(&iop, RBI(real)); + mpr_init(bigr); + mpr_set(bigr, bigratio); + + mpi_mul(&iop, &iop, mpr_den(bigratio)); + + if (nop == NOP_ADD) + mpi_add(mpr_num(bigr), &iop, mpr_num(bigr)); + else + mpi_sub(mpr_num(bigr), &iop, mpr_num(bigr)); + mpi_clear(&iop); + + RCLEAR_BI(real); + RBR(real) = bigr; + RTYPE(real) = N_BIGRATIO; + rbr_canonicalize(real); +} + +static void +rop_bi_br_md_xr(n_real *real, mpr *bigratio, int nop) +{ + mpr *bigr = XALLOC(mpr); + + mpr_init(bigr); + if (nop == NOP_MUL) + mpr_set(bigr, bigratio); + else + mpr_inv(bigr, bigratio); + + mpi_mul(mpr_num(bigr), RBI(real), mpr_num(bigr)); + + RCLEAR_BI(real); + RBR(real) = bigr; + RTYPE(real) = N_BIGRATIO; + rbr_canonicalize(real); +} + +static INLINE void +radd_bi_br(n_real *real, mpr *bigratio) +{ + rop_bi_br_as_xr(real, bigratio, NOP_ADD); +} + +static INLINE void +rsub_bi_br(n_real *real, mpr *bigratio) +{ + rop_bi_br_as_xr(real, bigratio, NOP_SUB); +} + +static INLINE void +rmul_bi_br(n_real *real, mpr *bigratio) +{ + rop_bi_br_md_xr(real, bigratio, NOP_MUL); +} + +static INLINE void +rdiv_bi_br(n_real *real, mpr *bigratio) +{ + rop_bi_br_md_xr(real, bigratio, NOP_DIV); +} + +static int +cmp_bi_br(mpi *bignum, mpr *bigratio) +{ + int cmp; + mpr cmp1; + + mpr_init(&cmp1); + mpi_set(mpr_num(&cmp1), bignum); + mpi_seti(mpr_den(&cmp1), 1); + + cmp = mpr_cmp(&cmp1, bigratio); + mpr_clear(&cmp1); + + return (cmp); +} + + +/************************************************************************ + * FIXRATIO FIXNUM + ************************************************************************/ +static void +rop_fr_fi_as_xr(n_real *real, long op, int nop) +{ + int fit; + long value = 0, num = RFRN(real), den = RFRD(real); + + fit = !fi_fi_mul_overflow(op, den); + + if (fit) { + value = op * den; + if (nop == NOP_ADD) + fit = !fi_fi_add_overflow(value, num); + else + fit = !fi_fi_sub_overflow(value, num); + } + if (fit) { + if (nop == NOP_ADD) + RFRN(real) = num + value; + else + RFRN(real) = num - value; + rfr_canonicalize(real); + } + else { + mpi iop; + mpr *bigr = XALLOC(mpr); + + mpr_init(bigr); + mpr_seti(bigr, num, den); + mpi_init(&iop); + mpi_seti(&iop, op); + mpi_muli(&iop, &iop, den); + if (nop == NOP_ADD) + mpi_add(mpr_num(bigr), mpr_num(bigr), &iop); + else + mpi_sub(mpr_num(bigr), mpr_num(bigr), &iop); + mpi_clear(&iop); + RBR(real) = bigr; + RTYPE(real) = N_BIGRATIO; + rbr_canonicalize(real); + } +} + +static void +rop_fr_fi_md_xr(n_real *real, long op, int nop) +{ + long num = RFRN(real), den = RFRD(real); + + if (nop == NOP_MUL) { + if (!fi_fi_mul_overflow(op, num)) { + RFRN(real) = op * num; + rfr_canonicalize(real); + return; + } + } + else if (!fi_fi_mul_overflow(op, den)) { + RFRD(real) = op * den; + rfr_canonicalize(real); + return; + } + + { + mpr *bigr = XALLOC(mpr); + + mpr_init(bigr); + mpr_seti(bigr, num, den); + if (nop == NOP_MUL) + mpr_muli(bigr, bigr, op); + else + mpr_divi(bigr, bigr, op); + RBR(real) = bigr; + RTYPE(real) = N_BIGRATIO; + rbr_canonicalize(real); + } +} + +static INLINE void +radd_fr_fi(n_real *real, long op) +{ + rop_fr_fi_as_xr(real, op, NOP_ADD); +} + +static INLINE void +rsub_fr_fi(n_real *real, long op) +{ + rop_fr_fi_as_xr(real, op, NOP_SUB); +} + +static INLINE void +rmul_fr_fi(n_real *real, long op) +{ + rop_fr_fi_md_xr(real, op, NOP_MUL); +} + +static INLINE void +rdiv_fr_fi(n_real *real, long op) +{ + rop_fr_fi_md_xr(real, op, NOP_DIV); +} + +static INLINE int +cmp_fr_fi(long num, long den, long fixnum) +{ + return (cmp_flonum((double)num / (double)den, (double)fixnum)); +} + + +/************************************************************************ + * FIXRATIO BIGNUM + ************************************************************************/ +static void +rop_fr_bi_as_xr(n_real *real, mpi *bignum, int nop) +{ + mpi iop; + mpr *bigr = XALLOC(mpr); + + mpr_init(bigr); + mpr_seti(bigr, RFRN(real), RFRD(real)); + + mpi_init(&iop); + mpi_set(&iop, bignum); + mpi_muli(&iop, &iop, RFRD(real)); + + if (nop == NOP_ADD) + mpi_add(mpr_num(bigr), mpr_num(bigr), &iop); + else + mpi_sub(mpr_num(bigr), mpr_num(bigr), &iop); + mpi_clear(&iop); + + RBR(real) = bigr; + RTYPE(real) = N_BIGRATIO; + rbr_canonicalize(real); +} + +static void +rop_fr_bi_md_xr(n_real *real, mpi *bignum, int nop) +{ + mpr *bigr = XALLOC(mpr); + + mpr_init(bigr); + mpr_seti(bigr, RFRN(real), RFRD(real)); + + if (nop == NOP_MUL) + mpi_mul(mpr_num(bigr), mpr_num(bigr), bignum); + else + mpi_mul(mpr_den(bigr), mpr_den(bigr), bignum); + + RBR(real) = bigr; + RTYPE(real) = N_BIGRATIO; + rbr_canonicalize(real); +} + +static INLINE void +radd_fr_bi(n_real *real, mpi *bignum) +{ + rop_fr_bi_as_xr(real, bignum, NOP_ADD); +} + +static INLINE void +rsub_fr_bi(n_real *real, mpi *bignum) +{ + rop_fr_bi_as_xr(real, bignum, NOP_SUB); +} + +static INLINE void +rmul_fr_bi(n_real *real, mpi *bignum) +{ + rop_fr_bi_md_xr(real, bignum, NOP_MUL); +} + +static INLINE void +rdiv_fr_bi(n_real *real, mpi *bignum) +{ + rop_fr_bi_md_xr(real, bignum, NOP_DIV); +} + +static int +cmp_fr_bi(long num, long den, mpi *bignum) +{ + int cmp; + mpr cmp1, cmp2; + + mpr_init(&cmp1); + mpr_seti(&cmp1, num, den); + + mpr_init(&cmp2); + mpi_set(mpr_num(&cmp2), bignum); + mpi_seti(mpr_den(&cmp2), 1); + + cmp = mpr_cmp(&cmp1, &cmp2); + mpr_clear(&cmp1); + mpr_clear(&cmp2); + + return (cmp); +} + + +/************************************************************************ + * FIXRATIO FIXRATIO + ************************************************************************/ +static void +rop_fr_fr_as_xr(n_real *real, long num2, long den2, int nop) +{ + int fit; + long num1 = RFRN(real), den1 = RFRD(real), num = 0, den = 0; + + fit = !fi_fi_mul_overflow(num1, den2); + if (fit) { + num = num1 * den2; + fit = !fi_fi_mul_overflow(num2, den1); + if (fit) { + den = num2 * den1; + if (nop == NOP_ADD) { + if ((fit = !fi_fi_add_overflow(num, den)) != 0) + num += den; + } + else if ((fit = !fi_fi_sub_overflow(num, den)) != 0) + num -= den; + if (fit) { + fit = !fi_fi_mul_overflow(den1, den2); + if (fit) + den = den1 * den2; + } + } + } + if (fit) { + RFRN(real) = num; + RFRD(real) = den; + rfr_canonicalize(real); + } + else { + mpi iop; + mpr *bigr = XALLOC(mpr); + + mpr_init(bigr); + mpr_seti(bigr, num1, den1); + mpi_muli(mpr_den(bigr), mpr_den(bigr), den2); + mpi_init(&iop); + mpi_seti(&iop, num2); + mpi_muli(&iop, &iop, den1); + mpi_muli(mpr_num(bigr), mpr_num(bigr), den2); + if (nop == NOP_ADD) + mpi_add(mpr_num(bigr), mpr_num(bigr), &iop); + else + mpi_sub(mpr_num(bigr), mpr_num(bigr), &iop); + mpi_clear(&iop); + RBR(real) = bigr; + RTYPE(real) = N_BIGRATIO; + rbr_canonicalize(real); + } +} + +static void +rop_fr_fr_md_xr(n_real *real, long num2, long den2, int nop) +{ + int fit; + long num1 = RFRN(real), den1 = RFRD(real), num = 0, den = 0; + + if (nop == NOP_MUL) { + fit = !fi_fi_mul_overflow(num1, num2) && !fi_fi_mul_overflow(den1, den2); + if (fit) { + num = num1 * num2; + den = den1 * den2; + } + } + else { + fit = !fi_fi_mul_overflow(num1, den2) && !fi_fi_mul_overflow(den1, num2); + if (fit) { + num = num1 * den2; + den = den1 * num2; + } + } + + if (fit) { + RFRN(real) = num; + RFRD(real) = den; + rfr_canonicalize(real); + } + else { + mpr *bigr = XALLOC(mpr); + + mpr_init(bigr); + + if (nop == NOP_MUL) { + mpr_seti(bigr, num1, den1); + mpi_muli(mpr_num(bigr), mpr_num(bigr), num2); + mpi_muli(mpr_den(bigr), mpr_den(bigr), den2); + } + else { + mpr_seti(bigr, num1, num2); + mpi_muli(mpr_num(bigr), mpr_num(bigr), den2); + mpi_muli(mpr_den(bigr), mpr_den(bigr), den1); + } + + RBR(real) = bigr; + RTYPE(real) = N_BIGRATIO; + rbr_canonicalize(real); + } +} + +static INLINE void +radd_fr_fr(n_real *real, long num, long den) +{ + rop_fr_fr_as_xr(real, num, den, NOP_ADD); +} + +static INLINE void +rsub_fr_fr(n_real *real, long num, long den) +{ + rop_fr_fr_as_xr(real, num, den, NOP_SUB); +} + +static INLINE void +rmul_fr_fr(n_real *real, long num, long den) +{ + rop_fr_fr_md_xr(real, num, den, NOP_MUL); +} + +static INLINE void +rdiv_fr_fr(n_real *real, long num, long den) +{ + rop_fr_fr_md_xr(real, num, den, NOP_DIV); +} + +static INLINE int +cmp_fr_fr(long num1, long den1, long num2, long den2) +{ + return (cmp_flonum((double)num1 / (double)den1, + (double)num2 / (double)den2)); +} + + +/************************************************************************ + * FIXRATIO BIGRATIO + ************************************************************************/ +static void +rop_fr_br_asmd_xr(n_real *real, mpr *bigratio, int nop) +{ + mpr *bigr = XALLOC(mpr); + + mpr_init(bigr); + mpr_seti(bigr, RFRN(real), RFRD(real)); + + switch (nop) { + case NOP_ADD: + mpr_add(bigr, bigr, bigratio); + break; + case NOP_SUB: + mpr_sub(bigr, bigr, bigratio); + break; + case NOP_MUL: + mpr_mul(bigr, bigr, bigratio); + break; + default: + mpr_div(bigr, bigr, bigratio); + break; + } + + RBR(real) = bigr; + RTYPE(real) = N_BIGRATIO; + rbr_canonicalize(real); +} + +static INLINE void +radd_fr_br(n_real *real, mpr *bigratio) +{ + rop_fr_br_asmd_xr(real, bigratio, NOP_ADD); +} + +static INLINE void +rsub_fr_br(n_real *real, mpr *bigratio) +{ + rop_fr_br_asmd_xr(real, bigratio, NOP_SUB); +} + +static INLINE void +rmul_fr_br(n_real *real, mpr *bigratio) +{ + rop_fr_br_asmd_xr(real, bigratio, NOP_MUL); +} + +static INLINE void +rdiv_fr_br(n_real *real, mpr *bigratio) +{ + rop_fr_br_asmd_xr(real, bigratio, NOP_DIV); +} + +static int +cmp_fr_br(long num, long den, mpr *bigratio) +{ + int cmp; + mpr cmp1; + + mpr_init(&cmp1); + mpr_seti(&cmp1, num, den); + + cmp = mpr_cmp(&cmp1, bigratio); + mpr_clear(&cmp1); + + return (cmp); +} + + +/************************************************************************ + * BIGRATIO FIXNUM + ************************************************************************/ +static void +rop_br_fi_asmd_xr(n_real *real, long fixnum, int nop) +{ + mpr *bigratio = RBR(real); + + switch (nop) { + case NOP_ADD: + mpr_addi(bigratio, bigratio, fixnum); + break; + case NOP_SUB: + mpr_subi(bigratio, bigratio, fixnum); + break; + case NOP_MUL: + mpr_muli(bigratio, bigratio, fixnum); + break; + default: + if (fixnum == 0) + fatal_error(DIVIDE_BY_ZERO); + mpr_divi(bigratio, bigratio, fixnum); + break; + } + rbr_canonicalize(real); +} + +static INLINE void +radd_br_fi(n_real *real, long fixnum) +{ + rop_br_fi_asmd_xr(real, fixnum, NOP_ADD); +} + +static INLINE void +rsub_br_fi(n_real *real, long fixnum) +{ + rop_br_fi_asmd_xr(real, fixnum, NOP_SUB); +} + +static INLINE void +rmul_br_fi(n_real *real, long fixnum) +{ + rop_br_fi_asmd_xr(real, fixnum, NOP_MUL); +} + +static INLINE void +rdiv_br_fi(n_real *real, long fixnum) +{ + rop_br_fi_asmd_xr(real, fixnum, NOP_DIV); +} + +static int +cmp_br_fi(mpr *bigratio, long fixnum) +{ + int cmp; + mpr cmp2; + + mpr_init(&cmp2); + mpr_seti(&cmp2, fixnum, 1); + cmp = mpr_cmp(bigratio, &cmp2); + mpr_clear(&cmp2); + + return (cmp); +} + + +/************************************************************************ + * BIGRATIO BIGNUM + ************************************************************************/ +static void +rop_br_bi_as_xr(n_real *real, mpi *bignum, int nop) +{ + mpi iop; + + mpi_init(&iop); + mpi_set(&iop, bignum); + + mpi_mul(&iop, &iop, RBRD(real)); + if (nop == NOP_ADD) + mpi_add(RBRN(real), RBRN(real), &iop); + else + mpi_sub(RBRN(real), RBRN(real), &iop); + mpi_clear(&iop); + rbr_canonicalize(real); +} + +static INLINE void +radd_br_bi(n_real *real, mpi *bignum) +{ + rop_br_bi_as_xr(real, bignum, NOP_ADD); +} + +static INLINE void +rsub_br_bi(n_real *real, mpi *bignum) +{ + rop_br_bi_as_xr(real, bignum, NOP_SUB); +} + +static INLINE void +rmul_br_bi(n_real *real, mpi *bignum) +{ + mpi_mul(RBRN(real), RBRN(real), bignum); + rbr_canonicalize(real); +} + +static INLINE void +rdiv_br_bi(n_real *real, mpi *bignum) +{ + mpi_mul(RBRD(real), RBRD(real), bignum); + rbr_canonicalize(real); +} + +static int +cmp_br_bi(mpr *bigratio, mpi *bignum) +{ + int cmp; + mpr cmp1; + + mpr_init(&cmp1); + mpi_set(mpr_num(&cmp1), bignum); + mpi_seti(mpr_den(&cmp1), 1); + + cmp = mpr_cmp(bigratio, &cmp1); + mpr_clear(&cmp1); + + return (cmp); +} + + +/************************************************************************ + * BIGRATIO FIXRATIO + ************************************************************************/ +static void +rop_br_fr_asmd_xr(n_real *real, long num, long den, int nop) +{ + mpr *bigratio = RBR(real), rop; + + mpr_init(&rop); + mpr_seti(&rop, num, den); + switch (nop) { + case NOP_ADD: + mpr_add(bigratio, bigratio, &rop); + break; + case NOP_SUB: + mpr_sub(bigratio, bigratio, &rop); + break; + case NOP_MUL: + mpr_mul(bigratio, bigratio, &rop); + break; + default: + mpr_div(bigratio, bigratio, &rop); + break; + } + mpr_clear(&rop); + rbr_canonicalize(real); +} + +static INLINE void +radd_br_fr(n_real *real, long num, long den) +{ + rop_br_fr_asmd_xr(real, num, den, NOP_ADD); +} + +static INLINE void +rsub_br_fr(n_real *real, long num, long den) +{ + rop_br_fr_asmd_xr(real, num, den, NOP_SUB); +} + +static INLINE void +rmul_br_fr(n_real *real, long num, long den) +{ + rop_br_fr_asmd_xr(real, num, den, NOP_MUL); +} + +static INLINE void +rdiv_br_fr(n_real *real, long num, long den) +{ + rop_br_fr_asmd_xr(real, num, den, NOP_DIV); +} + +static int +cmp_br_fr(mpr *bigratio, long num, long den) +{ + int cmp; + mpr cmp2; + + mpr_init(&cmp2); + mpr_seti(&cmp2, num, den); + cmp = mpr_cmp(bigratio, &cmp2); + mpr_clear(&cmp2); + + return (cmp); +} + + +/************************************************************************ + * BIGRATIO BIGRATIO + ************************************************************************/ +static INLINE void +radd_br_br(n_real *real, mpr *bigratio) +{ + mpr_add(RBR(real), RBR(real), bigratio); + rbr_canonicalize(real); +} + +static INLINE void +rsub_br_br(n_real *real, mpr *bigratio) +{ + mpr_sub(RBR(real), RBR(real), bigratio); + rbr_canonicalize(real); +} + +static INLINE void +rmul_br_br(n_real *real, mpr *bigratio) +{ + mpr_mul(RBR(real), RBR(real), bigratio); + rbr_canonicalize(real); +} + +static INLINE void +rdiv_br_br(n_real *real, mpr *bigratio) +{ + mpr_div(RBR(real), RBR(real), bigratio); + rbr_canonicalize(real); +} + +static INLINE int +cmp_br_br(mpr *op1, mpr *op2) +{ + return (mpr_cmp(op1, op2)); +} |