6644N/AA new version of mpfr is expected which isn't in Solaris yet, so we
6644N/Aare removing that functionality till we can get mpfr updated and
6644N/Athen will revert back these changes.
6644N/A
6644N/ANot suitable for upstream
6644N/A
6644N/A--- a/data/org.gnome.calculator.gschema.xml Fri Jan 29 15:39:20 2016
6644N/A+++ b/data/org.gnome.calculator.gschema.xml Fri Jan 29 15:40:20 2016
6644N/A@@ -22,6 +22,7 @@
6644N/A <schema path="/org/gnome/calculator/" id="org.gnome.calculator" gettext-domain="calculator">
6644N/A <key type="i" name="accuracy">
6644N/A <default>9</default>
6644N/A+ <range min="0" max="9"/>
6644N/A <summary>Accuracy value</summary>
6644N/A <description>The number of digits displayed after the numeric point</description>
6644N/A </key>
6644N/A@@ -82,10 +83,5 @@
6644N/A <summary>Target units</summary>
6644N/A <description>Units to convert the current calculation into</description>
6644N/A </key>
6644N/A- <key type="i" name="precision">
6644N/A- <default>2000</default>
6644N/A- <summary>Internal precision</summary>
6644N/A- <description>The internal precision used with the MPFR library</description>
6644N/A- </key>
6644N/A </schema>
6644N/A </schemalist>
6644N/A--- a/src/Makefile.am Fri Jan 29 15:40:46 2016
6644N/A+++ b/src/Makefile.am Fri Jan 29 15:41:34 2016
6644N/A@@ -33,8 +33,7 @@
6644N/A --pkg gtk+-3.0 \
6644N/A --pkg gtksourceview-3.0 \
6644N/A --pkg libxml-2.0 \
6644N/A- $(top_builddir)/lib/libcalculator.vapi \
6644N/A- $(top_builddir)/lib/mpfr.vapi
6644N/A+ $(top_builddir)/lib/libcalculator.vapi
6644N/A
6644N/A gnome_calculator_CPPFLAGS = \
6644N/A $(AM_CPPFLAGS) \
6644N/A@@ -54,8 +53,7 @@
6644N/A --pkg gio-2.0 \
6644N/A --pkg gtksourceview-3.0 \
6644N/A --pkg libxml-2.0 \
6644N/A- $(top_builddir)/lib/libcalculator.vapi \
6644N/A- $(top_builddir)/lib/mpfr.vapi
6644N/A+ $(top_builddir)/lib/libcalculator.vapi
6644N/A
6644N/A gcalccmd_LDADD = \
6644N/A $(top_builddir)/lib/libcalculator.la
6644N/A--- a/lib/Makefile.am Fri Jan 29 15:42:13 2016
6644N/A+++ b/lib/Makefile.am Fri Jan 29 15:42:24 2016
6644N/A@@ -10,7 +10,6 @@
6644N/A
6644N/A libcalculator_la_SOURCES = \
6644N/A config.vapi \
6644N/A- mpfr.vapi \
6644N/A currency.vala \
6644N/A equation.vala \
6644N/A equation-lexer.vala \
6644N/A@@ -37,8 +36,7 @@
6644N/A libcalculator_la_LIBADD = \
6644N/A $(LIBCALCULATOR_LIBS) \
6644N/A -lgmp \
6644N/A- -lm \
6644N/A- -lmpfr
6644N/A+ -lm
6644N/A
6644N/A EXTRA_DIST = \
6644N/A libcalculator.h \
6644N/A--- a/lib/equation-parser.vala Fri Jan 29 15:47:35 2016
6644N/A+++ b/lib/equation-parser.vala Fri Jan 29 16:00:51 2016
6644N/A@@ -78,20 +78,8 @@
6644N/A var r = right.solve ();
6644N/A if (r == null)
6644N/A return null;
6644N/A- var z = solve_r (r);
6644N/A+ return solve_r (r);
6644N/A
6644N/A- /* check for errors */
6644N/A- Number.check_flags ();
6644N/A- if (Number.error != null)
6644N/A- {
6644N/A- var tmpleft = right;
6644N/A- var tmpright = right;
6644N/A- while (tmpleft.left != null) tmpleft = tmpleft.left;
6644N/A- while (tmpright.right != null) tmpright = tmpright.right;
6644N/A- parser.set_error (ErrorCode.MP, Number.error, tmpleft.token.start_index, tmpright.token.end_index);
6644N/A- Number.error = null;
6644N/A- }
6644N/A- return z;
6644N/A }
6644N/A
6644N/A public abstract Number solve_r (Number r);
6644N/A@@ -110,20 +98,7 @@
6644N/A var r = right.solve ();
6644N/A if (l == null || r == null)
6644N/A return null;
6644N/A- var z = solve_lr (l, r);
6644N/A-
6644N/A- /* check for errors */
6644N/A- Number.check_flags ();
6644N/A- if (Number.error != null)
6644N/A- {
6644N/A- var tmpleft = left;
6644N/A- var tmpright = right;
6644N/A- while (tmpleft.left != null) tmpleft = tmpleft.left;
6644N/A- while (tmpright.right != null) tmpright = tmpright.right;
6644N/A- parser.set_error (ErrorCode.MP, Number.error, tmpleft.token.start_index, tmpright.token.end_index);
6644N/A- Number.error = null;
6644N/A- }
6644N/A- return z;
6644N/A+ return solve_lr (l, r);
6644N/A }
6644N/A
6644N/A public abstract Number solve_lr (Number left, Number r);
6644N/A@@ -261,18 +236,6 @@
6644N/A value = value.multiply (t);
6644N/A }
6644N/A
6644N/A- /* check for errors */
6644N/A- Number.check_flags ();
6644N/A- if (Number.error != null)
6644N/A- {
6644N/A- var tmpleft = left;
6644N/A- var tmpright = right;
6644N/A- while (tmpleft.left != null) tmpleft = tmpleft.left;
6644N/A- while (tmpright.right != null) tmpright = tmpright.right;
6644N/A- parser.set_error (ErrorCode.MP, Number.error, tmpleft.token.start_index, tmpright.token.end_index);
6644N/A- Number.error = null;
6644N/A- }
6644N/A-
6644N/A return value;
6644N/A }
6644N/A }
6644N/A@@ -392,14 +355,6 @@
6644N/A if (tmp != null)
6644N/A tmp = tmp.xpowy_integer (pow);
6644N/A
6644N/A- /* check for errors */
6644N/A- Number.check_flags ();
6644N/A- if (Number.error != null)
6644N/A- {
6644N/A- parser.set_error (ErrorCode.MP, Number.error);
6644N/A- Number.error = null;
6644N/A- }
6644N/A-
6644N/A return tmp;
6644N/A }
6644N/A }
6644N/A@@ -574,17 +529,7 @@
6644N/A
6644N/A public override Number solve_lr (Number l, Number r)
6644N/A {
6644N/A- var z = l.divide (r);
6644N/A- if (Number.error != null)
6644N/A- {
6644N/A- var tmpleft = right;
6644N/A- var tmpright = right;
6644N/A- while (tmpleft.left != null) tmpleft = tmpleft.left;
6644N/A- while (tmpright.right != null) tmpright = tmpright.right;
6644N/A- parser.set_error (ErrorCode.MP, Number.error, tmpleft.token.start_index, tmpright.token.end_index);
6644N/A- Number.error = null;
6644N/A- }
6644N/A- return z;
6644N/A+ return l.divide (r);
6644N/A }
6644N/A }
6644N/A
6644N/A@@ -604,21 +549,8 @@
6644N/A var mod = right.solve ();
6644N/A if (base_value == null || exponent == null || mod == null)
6644N/A return null;
6644N/A- var z = base_value.modular_exponentiation (exponent, mod);
6644N/A+ return base_value.modular_exponentiation (exponent, mod);
6644N/A
6644N/A- /* check for errors */
6644N/A- Number.check_flags ();
6644N/A- if (Number.error != null)
6644N/A- {
6644N/A- var tmpleft = left;
6644N/A- var tmpright = right;
6644N/A- while (tmpleft.left != null) tmpleft = tmpleft.left;
6644N/A- while (tmpright.right != null) tmpright = tmpright.right;
6644N/A- parser.set_error (ErrorCode.MP, Number.error, tmpleft.token.start_index, tmpright.token.end_index);
6644N/A- Number.error = null;
6644N/A- }
6644N/A-
6644N/A- return z;
6644N/A }
6644N/A else
6644N/A {
6644N/A@@ -626,21 +558,8 @@
6644N/A var r = right.solve ();
6644N/A if (l == null || r == null)
6644N/A return null;
6644N/A- var z = solve_lr (l, r);
6644N/A+ return solve_lr (l, r);
6644N/A
6644N/A- /* check for errors */
6644N/A- Number.check_flags ();
6644N/A- if (Number.error != null)
6644N/A- {
6644N/A- var tmpleft = left;
6644N/A- var tmpright = right;
6644N/A- while (tmpleft.left != null) tmpleft = tmpleft.left;
6644N/A- while (tmpright.right != null) tmpright = tmpright.right;
6644N/A- parser.set_error (ErrorCode.MP, Number.error, tmpleft.token.start_index, tmpright.token.end_index);
6644N/A- Number.error = null;
6644N/A- }
6644N/A-
6644N/A- return z;
6644N/A }
6644N/A }
6644N/A
6644N/A@@ -709,21 +628,8 @@
6644N/A if (val == null)
6644N/A return null;
6644N/A
6644N/A- var z = val.xpowy_integer (pow);
6644N/A+ return val.xpowy_integer (pow);
6644N/A
6644N/A- /* check for errors */
6644N/A- Number.check_flags ();
6644N/A- if (Number.error != null)
6644N/A- {
6644N/A- var tmpleft = left;
6644N/A- var tmpright = right;
6644N/A- while (tmpleft.left != null) tmpleft = tmpleft.left;
6644N/A- while (tmpright.right != null) tmpright = tmpright.right;
6644N/A- parser.set_error (ErrorCode.MP, Number.error, tmpleft.token.start_index, tmpright.token.end_index);
6644N/A- Number.error = null;
6644N/A- }
6644N/A-
6644N/A- return z;
6644N/A }
6644N/A }
6644N/A
6644N/A@@ -1036,10 +942,10 @@
6644N/A }
6644N/A
6644N/A representation_base = this.representation_base;
6644N/A- error_code = this.error;
6644N/A- error_token = this.error_token;
6644N/A- error_start = this.error_token_start;
6644N/A- error_end = this.error_token_end;
6644N/A+ error_code = ErrorCode.NONE;
6644N/A+ error_token = null;
6644N/A+ error_start = 0;
6644N/A+ error_end = 0;
6644N/A return ans;
6644N/A }
6644N/A
6644N/A--- a/lib/equation.vala Fri Jan 29 16:01:11 2016
6644N/A+++ b/lib/equation.vala Fri Jan 29 16:02:09 2016
6644N/A@@ -117,17 +117,15 @@
6644N/A public new Number? parse (out uint representation_base = null, out ErrorCode error_code = null, out string? error_token = null, out uint? error_start = null, out uint? error_end = null)
6644N/A {
6644N/A var parser = new EquationParser (this, expression);
6644N/A- Number.error = null;
6644N/A+ mp_clear_error();
6644N/A
6644N/A var z = parser.parse (out representation_base, out error_code, out error_token, out error_start, out error_end);
6644N/A
6644N/A /* Error during parsing */
6644N/A if (error_code != ErrorCode.NONE)
6644N/A- {
6644N/A return null;
6644N/A- }
6644N/A
6644N/A- if (Number.error != null)
6644N/A+ if (mp_get_error() != null)
6644N/A {
6644N/A error_code = ErrorCode.MP;
6644N/A return null;
6644N/A--- a/src/gnome-calculator.vala Fri Jan 29 16:05:17 2016
6644N/A+++ b/src/gnome-calculator.vala Fri Jan 29 16:06:27 2016
6644N/A@@ -56,7 +56,6 @@
6644N/A var target_currency = settings.get_string ("target-currency");
6644N/A var source_units = settings.get_string ("source-units");
6644N/A var target_units = settings.get_string ("target-units");
6644N/A- var precision = settings.get_int ("precision");
6644N/A
6644N/A var equation = new MathEquation ();
6644N/A equation.accuracy = accuracy;
6644N/A@@ -69,7 +68,6 @@
6644N/A equation.target_currency = target_currency;
6644N/A equation.source_units = source_units;
6644N/A equation.target_units = target_units;
6644N/A- Number.precision = precision;
6644N/A
6644N/A add_action_entries (app_entries, this);
6644N/A
6644N/A@@ -173,7 +171,7 @@
6644N/A }
6644N/A else if (error == ErrorCode.MP)
6644N/A {
6644N/A- stderr.printf ("Error: %s\n", Number.error);
6644N/A+ stderr.printf ("Error: %s\n", mp_get_error());
6644N/A return Posix.EXIT_FAILURE;
6644N/A }
6644N/A else
6644N/A--- a/src/gcalccmd.vala Fri Jan 29 16:03:42 2016
6644N/A+++ b/src/gcalccmd.vala Fri Jan 29 16:04:56 2016
6644N/A@@ -34,18 +34,9 @@
6644N/A
6644N/A result_serializer.set_representation_base (representation_base);
6644N/A if (z != null)
6644N/A- {
6644N/A- var str = result_serializer.to_string (z);
6644N/A- if (result_serializer.error != null)
6644N/A- {
6644N/A- stderr.printf ("%s\n", result_serializer.error);
6644N/A- result_serializer.error = null;
6644N/A- }
6644N/A- else
6644N/A- stdout.printf ("%s\n", str);
6644N/A- }
6644N/A+ stdout.printf ("%s\n", result_serializer.to_string (z));
6644N/A else if (ret == ErrorCode.MP)
6644N/A- stderr.printf ("Error %s\n", Number.error);
6644N/A+ stderr.printf ("Error %s\n", mp_get_error());
6644N/A else
6644N/A stderr.printf ("Error %d\n", ret);
6644N/A }
6644N/A--- a/lib/math-equation.vala Fri Jan 29 16:33:36 2016
6644N/A+++ b/lib/math-equation.vala Fri Jan 29 16:39:56 2016
6644N/A@@ -786,11 +786,6 @@
6644N/A ans_end_mark = create_mark (null, end, true);
6644N/A apply_tag (ans_tag, start, end);
6644N/A
6644N/A- if (serializer.error != null)
6644N/A- {
6644N/A- status = serializer.error;
6644N/A- serializer.error = null;
6644N/A- }
6644N/A }
6644N/A
6644N/A public new void insert (string text)
6644N/A@@ -964,13 +959,11 @@
6644N/A break;
6644N/A
6644N/A case ErrorCode.MP:
6644N/A- if (Number.error != null) // LEGACY, should never be run
6644N/A+ if (mp_get_error () != null)
6644N/A+ solvedata.error = mp_get_error ();
6644N/A+ else if (error_token != null)
6644N/A {
6644N/A- solvedata.error = Number.error;
6644N/A- }
6644N/A- else if (error_token != null) // should always be run
6644N/A- {
6644N/A- solvedata.error = _("%s").printf (error_token);
6644N/A+ solvedata.error = _("Malformed expression at token '%s'").printf(error_token);
6644N/A solvedata.error_start = error_start;
6644N/A solvedata.error_end = error_end;
6644N/A }
6644N/A@@ -1015,8 +1008,6 @@
6644N/A
6644N/A /* Fix thousand separator offsets in the start and end offsets of error token. */
6644N/A error_token_fix_thousands_separator ();
6644N/A- /* Fix missing Parenthesis before the start and after the end offsets of error token */
6644N/A- error_token_fix_parenthesis ();
6644N/A
6644N/A /* Notify the GUI about the change in error token locations. */
6644N/A notify_property ("error-token-end");
6644N/A@@ -1087,70 +1078,6 @@
6644N/A end.forward_chars (length);
6644N/A }
6644N/A }
6644N/A-
6644N/A- /* Fix the offsets to consider starting and ending parenthesis */
6644N/A- private void error_token_fix_parenthesis ()
6644N/A- {
6644N/A- unichar c;
6644N/A- int count = 0;
6644N/A- int real_end = display.index_of_nth_char (error_token_end);
6644N/A- int real_start = display.index_of_nth_char (error_token_start);
6644N/A-
6644N/A- /* checks if there are more opening/closing parenthesis than closing/opening parenthesis */
6644N/A- for (int i = real_start; display.get_next_char (ref i, out c) && i <= real_end;)
6644N/A- {
6644N/A- if (c.to_string () == "(") count++;
6644N/A- if (c.to_string () == ")") count--;
6644N/A- }
6644N/A-
6644N/A- /* if there are more opening than closing parenthesis and there are closing parenthesis
6644N/A- after the end offset, include those in the offsets */
6644N/A- for (int i = real_end; display.get_next_char (ref i, out c) && count > 0;)
6644N/A- {
6644N/A- if (c.to_string () == ")")
6644N/A- {
6644N/A- state.error_token_end++;
6644N/A- count--;
6644N/A- }
6644N/A- else
6644N/A- {
6644N/A- break;
6644N/A- }
6644N/A- }
6644N/A-
6644N/A- /* the same for closing parenthesis */
6644N/A- for (int i = real_start; display.get_prev_char (ref i, out c) && count < 0;)
6644N/A- {
6644N/A- if (c.to_string () == "(")
6644N/A- {
6644N/A- state.error_token_start--;
6644N/A- count++;
6644N/A- }
6644N/A- else
6644N/A- {
6644N/A- break;
6644N/A- }
6644N/A- }
6644N/A-
6644N/A- real_end = display.index_of_nth_char (error_token_end);
6644N/A- real_start = display.index_of_nth_char (error_token_start);
6644N/A-
6644N/A- unichar d;
6644N/A-
6644N/A- /* if there are opening parenthesis directly before aswell as closing parenthesis directly after the offsets, include those aswell */
6644N/A- while (display.get_next_char (ref real_end, out d) && display.get_prev_char (ref real_start, out c))
6644N/A- {
6644N/A- if (c.to_string () == "(" && d.to_string () == ")")
6644N/A- {
6644N/A- state.error_token_start--;
6644N/A- state.error_token_end++;
6644N/A- }
6644N/A- else
6644N/A- {
6644N/A- break;
6644N/A- }
6644N/A- }
6644N/A- }
6644N/A
6644N/A private void* factorize_real ()
6644N/A {
6644N/A--- a/src/math-preferences.vala Fri Jan 29 16:40:37 2016
6644N/A+++ b/src/math-preferences.vala Fri Jan 29 16:43:30 2016
6644N/A@@ -90,7 +90,7 @@
6644N/A var string = _("Show %d decimal _places");
6644N/A var tokens = string.split ("%d", 2);
6644N/A
6644N/A- var decimal_places_adjustment = new Gtk.Adjustment (0.0, 0.0, 100.0, 1.0, 1.0, 0.0);
6644N/A+ var decimal_places_adjustment = new Gtk.Adjustment (0.0, 0.0, 9.0, 1.0, 1.0, 0.0);
6644N/A decimal_places_spin = new Gtk.SpinButton (decimal_places_adjustment, 0.0, 0);
6644N/A
6644N/A if (tokens.length > 0)
6644N/A--- a/lib/serializer.vala Mon Feb 1 10:45:23 2016
6644N/A+++ b/lib/serializer.vala Mon Feb 1 10:46:37 2016
6644N/A@@ -32,9 +32,6 @@
6644N/A private unichar tsep; /* Locale specific thousands separator. */
6644N/A private int tsep_count; /* Number of digits between separator. */
6644N/A
6644N/A- /* is set when an error (for example precision error while converting) occurs */
6644N/A- public string? error { get; set; default = null; }
6644N/A-
6644N/A public Serializer (DisplayFormat format, int number_base, int trailing_digits)
6644N/A {
6644N/A var radix_string = Posix.nl_langinfo (Posix.NLItem.RADIXCHAR);
6644N/A@@ -322,17 +319,8 @@
6644N/A
6644N/A var d = t3.to_integer ();
6644N/A
6644N/A- if (d < 16 && d >= 0)
6644N/A- {
6644N/A- string.prepend_c (digits[d]);
6644N/A- }
6644N/A- else
6644N/A- {
6644N/A- string.prepend_c ('?');
6644N/A- error = _("Precision error");
6644N/A- string.assign ("0");
6644N/A- break;
6644N/A- }
6644N/A+ string.prepend_c (d < 16 ? digits[d] : '?');
6644N/A+
6644N/A n_digits++;
6644N/A
6644N/A temp = t;
6644N/A--- a/tests/test-number.vala Mon Feb 1 10:48:08 2016
6644N/A+++ b/tests/test-number.vala Mon Feb 1 10:49:35 2016
6644N/A@@ -72,6 +72,22 @@
6644N/A pass ();
6644N/A }
6644N/A
6644N/A+private void test_float ()
6644N/A+{
6644N/A+ for (var a = -10.0f; a <= 10.0f; a += 0.5f)
6644N/A+ {
6644N/A+ var z = new Number.float (a);
6644N/A+ if (z.to_float () != a)
6644N/A+ {
6644N/A+ fail ("Number.float (%f).to_float () -> %f, expected %f".printf (a, z.to_float (), a));
6644N/A+ return;
6644N/A+ }
6644N/A+ }
6644N/A+
6644N/A+ pass ();
6644N/A+}
6644N/A+
6644N/A+
6644N/A private void test_double ()
6644N/A {
6644N/A for (var a = -10.0; a <= 10.0; a += 0.5)
6644N/A@@ -1074,6 +1090,7 @@
6644N/A test_integer ();
6644N/A test_unsigned_integer ();
6644N/A test_fraction ();
6644N/A+ test_float ();
6644N/A test_double ();
6644N/A test_complex ();
6644N/A test_polar ();
6644N/A--- a/tests/test-equation.c Mon Feb 1 12:22:33 2016
6644N/A+++ b/tests/test-equation.c Mon Feb 1 12:22:47 2016
6644N/A@@ -95,7 +95,7 @@
6644N/A const gchar* _tmp1_ = NULL;
6644N/A const gchar* _tmp2_ = NULL;
6644N/A gchar* _tmp3_ = NULL;
6644N/A- _tmp1_ = number_get_error ();
6644N/A+ _tmp1_ = mp_get_error ();
6644N/A _tmp2_ = _tmp1_;
6644N/A _tmp3_ = g_strdup_printf ("ErrorCode.MP(\"%s\")", _tmp2_);
6644N/A result = _tmp3_;
6644N/A
6644N/A--- a/lib/number.vala Fri Jan 29 16:44:36 2016
6644N/A+++ b/lib/number.vala Mon Feb 1 12:02:40 2016
6644N/A@@ -27,8 +27,15 @@
6644N/A * THE MP USERS GUIDE.
6644N/A */
6644N/A
6644N/A-using MPFR;
6644N/A
6644N/A+/* Size of the multiple precision values */
6644N/A+private const int SIZE = 1000;
6644N/A+
6644N/A+/* Base for numbers */
6644N/A+private const int BASE = 10000;
6644N/A+
6644N/A+private const int T = 100;
6644N/A+
6644N/A private delegate int BitwiseFunc (int v1, int v2);
6644N/A
6644N/A public enum AngleUnit
6644N/A@@ -41,37 +48,69 @@
6644N/A /* Object for a high precision floating point number representation */
6644N/A public class Number : Object
6644N/A {
6644N/A- /* real and imaginary part of a Number */
6644N/A- private MPFloat re_num { get; set; }
6644N/A- private MPFloat im_num { get; set; }
6644N/A+ /* Sign (+1, -1) or 0 for the value zero */
6644N/A+ public int re_sign;
6644N/A+ public int im_sign;
6644N/A
6644N/A- public static ulong precision { get; set; default = 1000; }
6644N/A+ /* Exponent (to base BASE) */
6644N/A+ public int re_exponent;
6644N/A+ public int im_exponent;
6644N/A
6644N/A- /* Stores the error msg if an error occurs during calculation. Otherwise should be null */
6644N/A- public static string? error { get; set; default = null; }
6644N/A+ /* Normalized fraction */
6644N/A+ public int re_fraction[1000]; // SIZE
6644N/A+ public int im_fraction[1000]; // SIZE
6644N/A
6644N/A public Number.integer (int64 value)
6644N/A {
6644N/A- var tmp = MPFloat.init2 ((Precision) precision);
6644N/A- tmp.set_signed_integer ((long)value, Round.NEAREST);
6644N/A- re_num = tmp;
6644N/A- var tmp2 = MPFloat.init2 ((Precision) precision);
6644N/A- tmp2.set_unsigned_integer (0, Round.NEAREST);
6644N/A- im_num = tmp2;
6644N/A+ if (value < 0)
6644N/A+ {
6644N/A+ value = -value;
6644N/A+ re_sign = -1;
6644N/A+ }
6644N/A+ else if (value > 0)
6644N/A+ re_sign = 1;
6644N/A+ else
6644N/A+ re_sign = 0;
6644N/A+
6644N/A+ while (value != 0)
6644N/A+ {
6644N/A+ re_fraction[re_exponent] = (int) (value % BASE);
6644N/A+ re_exponent++;
6644N/A+ value /= BASE;
6644N/A+ }
6644N/A+ for (var i = 0; i < re_exponent / 2; i++)
6644N/A+ {
6644N/A+ int t = re_fraction[i];
6644N/A+ re_fraction[i] = re_fraction[re_exponent - i - 1];
6644N/A+ re_fraction[re_exponent - i - 1] = t;
6644N/A+ }
6644N/A }
6644N/A
6644N/A public Number.unsigned_integer (uint64 x)
6644N/A {
6644N/A- var tmp = MPFloat.init2 ((Precision) precision);
6644N/A- tmp.set_unsigned_integer ((ulong)x, Round.NEAREST);
6644N/A- re_num = tmp;
6644N/A- var tmp2 = MPFloat.init2 ((Precision) precision);
6644N/A- tmp2.set_unsigned_integer (0, Round.NEAREST);
6644N/A- im_num = tmp2;
6644N/A+ if (x == 0)
6644N/A+ re_sign = 0;
6644N/A+ else
6644N/A+ re_sign = 1;
6644N/A+
6644N/A+ while (x != 0)
6644N/A+ {
6644N/A+ re_fraction[re_exponent] = (int) (x % BASE);
6644N/A+ x = x / BASE;
6644N/A+ re_exponent++;
6644N/A+ }
6644N/A+ for (var i = 0; i < re_exponent / 2; i++)
6644N/A+ {
6644N/A+ int t = re_fraction[i];
6644N/A+ re_fraction[i] = re_fraction[re_exponent - i - 1];
6644N/A+ re_fraction[re_exponent - i - 1] = t;
6644N/A+ }
6644N/A }
6644N/A
6644N/A public Number.fraction (int64 numerator, int64 denominator)
6644N/A {
6644N/A+ mp_gcd (ref numerator, ref denominator);
6644N/A+
6644N/A if (denominator < 0)
6644N/A {
6644N/A numerator = -numerator;
6644N/A@@ -81,39 +120,203 @@
6644N/A Number.integer (numerator);
6644N/A if (denominator != 1)
6644N/A {
6644N/A- var tmp = re_num;
6644N/A- tmp.divide_signed_integer (re_num, (long) denominator, Round.NEAREST);
6644N/A- re_num = tmp;
6644N/A+ var z = divide_integer (denominator);
6644N/A+ re_sign = z.re_sign;
6644N/A+ im_sign = z.im_sign;
6644N/A+ re_exponent = z.re_exponent;
6644N/A+ im_exponent = z.im_exponent;
6644N/A+ for (var i = 0; i < z.re_fraction.length; i++)
6644N/A+ {
6644N/A+ re_fraction[i] = z.re_fraction[i];
6644N/A+ im_fraction[i] = z.im_fraction[i];
6644N/A+ }
6644N/A }
6644N/A }
6644N/A
6644N/A- /* Helper constructor. Creates new Number from already existing MPFloat. */
6644N/A- public Number.mpfloat (MPFloat value)
6644N/A+ public Number.float (float value)
6644N/A {
6644N/A- re_num = value;
6644N/A- var tmp = MPFloat.init2 ((Precision) precision);
6644N/A- tmp.set_unsigned_integer (0, Round.NEAREST);
6644N/A- im_num = tmp;
6644N/A+ var z = new Number.integer (0);
6644N/A+
6644N/A+ if (value != 0)
6644N/A+ {
6644N/A+ /* CHECK SIGN */
6644N/A+ var rj = 0f;
6644N/A+ if (value < 0.0f)
6644N/A+ {
6644N/A+ z.re_sign = -1;
6644N/A+ rj = -value;
6644N/A+ }
6644N/A+ else if (value > 0.0f)
6644N/A+ {
6644N/A+ z.re_sign = 1;
6644N/A+ rj = value;
6644N/A+ }
6644N/A+
6644N/A+ /* INCREASE IE AND DIVIDE RJ BY 16. */
6644N/A+ var ie = 0;
6644N/A+ while (rj >= 1.0f)
6644N/A+ {
6644N/A+ ie++;
6644N/A+ rj *= 0.0625f;
6644N/A+ }
6644N/A+ while (rj < 0.0625f)
6644N/A+ {
6644N/A+ ie--;
6644N/A+ rj *= 16.0f;
6644N/A+ }
6644N/A+
6644N/A+ /* NOW RJ IS DY DIVIDED BY SUITABLE POWER OF 16.
6644N/A+ * SET re_exponent TO 0
6644N/A+ */
6644N/A+ z.re_exponent = 0;
6644N/A+
6644N/A+ /* CONVERSION LOOP (ASSUME SINGLE-PRECISION OPS. EXACT) */
6644N/A+ for (var i = 0; i < T + 4; i++)
6644N/A+ {
6644N/A+ rj *= BASE;
6644N/A+ z.re_fraction[i] = (int) rj;
6644N/A+ rj -= z.re_fraction[i];
6644N/A+ }
6644N/A+
6644N/A+ /* NORMALIZE RESULT */
6644N/A+ mp_normalize (ref z);
6644N/A+
6644N/A+ /* Computing MAX */
6644N/A+ var ib = int.max (BASE * 7 * BASE, 32767) / 16;
6644N/A+ var tp = 1;
6644N/A+
6644N/A+ /* NOW MULTIPLY BY 16**IE */
6644N/A+ if (ie < 0)
6644N/A+ {
6644N/A+ var k = -ie;
6644N/A+ for (var i = 1; i <= k; i++)
6644N/A+ {
6644N/A+ tp <<= 4;
6644N/A+ if (tp <= ib && tp != BASE && i < k)
6644N/A+ continue;
6644N/A+ z = z.divide_integer (tp);
6644N/A+ tp = 1;
6644N/A+ }
6644N/A+ }
6644N/A+ else if (ie > 0)
6644N/A+ {
6644N/A+ for (var i = 1; i <= ie; i++)
6644N/A+ {
6644N/A+ tp <<= 4;
6644N/A+ if (tp <= ib && tp != BASE && i < ie)
6644N/A+ continue;
6644N/A+ z = z.multiply_integer (tp);
6644N/A+ tp = 1;
6644N/A+ }
6644N/A+ }
6644N/A+ }
6644N/A+
6644N/A+ re_sign = z.re_sign;
6644N/A+ im_sign = z.im_sign;
6644N/A+ re_exponent = z.re_exponent;
6644N/A+ im_exponent = z.im_exponent;
6644N/A+ for (var i = 0; i < z.re_fraction.length; i++)
6644N/A+ {
6644N/A+ re_fraction[i] = z.re_fraction[i];
6644N/A+ im_fraction[i] = z.im_fraction[i];
6644N/A+ }
6644N/A }
6644N/A
6644N/A public Number.double (double value)
6644N/A {
6644N/A- var tmp = MPFloat.init2 ((Precision) precision);
6644N/A- tmp.set_double (value, Round.NEAREST);
6644N/A- re_num = tmp;
6644N/A- var tmp2 = MPFloat.init2 ((Precision) precision);
6644N/A- tmp2.set_unsigned_integer (0, Round.NEAREST);
6644N/A- im_num = tmp2;
6644N/A+ var z = new Number.integer (0);
6644N/A+
6644N/A+ if (value != 0)
6644N/A+ {
6644N/A+ /* CHECK SIGN */
6644N/A+ var dj = 0.0;
6644N/A+ if (value < 0.0)
6644N/A+ {
6644N/A+ z.re_sign = -1;
6644N/A+ dj = -value;
6644N/A+ }
6644N/A+ else if (value > 0.0)
6644N/A+ {
6644N/A+ z.re_sign = 1;
6644N/A+ dj = value;
6644N/A+ }
6644N/A+
6644N/A+ /* INCREASE IE AND DIVIDE DJ BY 16. */
6644N/A+ var ie = 0;
6644N/A+ for (ie = 0; dj >= 1.0; ie++)
6644N/A+ dj *= 1.0/16.0;
6644N/A+
6644N/A+ for ( ; dj < 1.0/16.0; ie--)
6644N/A+ dj *= 16.0;
6644N/A+
6644N/A+ /* NOW DJ IS DY DIVIDED BY SUITABLE POWER OF 16
6644N/A+ * SET re_exponent TO 0
6644N/A+ */
6644N/A+ z.re_exponent = 0;
6644N/A+
6644N/A+ /* CONVERSION LOOP (ASSUME DOUBLE-PRECISION OPS. EXACT) */
6644N/A+ for (var i = 0; i < T + 4; i++)
6644N/A+ {
6644N/A+ dj *= (double) BASE;
6644N/A+ z.re_fraction[i] = (int) dj;
6644N/A+ dj -= (double) z.re_fraction[i];
6644N/A+ }
6644N/A+
6644N/A+ /* NORMALIZE RESULT */
6644N/A+ mp_normalize (ref z);
6644N/A+
6644N/A+ /* Computing MAX */
6644N/A+ var ib = int.max (BASE * 7 * BASE, 32767) / 16;
6644N/A+ var tp = 1;
6644N/A+
6644N/A+ /* NOW MULTIPLY BY 16**IE */
6644N/A+ if (ie < 0)
6644N/A+ {
6644N/A+ var k = -ie;
6644N/A+ for (var i = 1; i <= k; ++i)
6644N/A+ {
6644N/A+ tp <<= 4;
6644N/A+ if (tp <= ib && tp != BASE && i < k)
6644N/A+ continue;
6644N/A+ z = z.divide_integer (tp);
6644N/A+ tp = 1;
6644N/A+ }
6644N/A+ }
6644N/A+ else if (ie > 0)
6644N/A+ {
6644N/A+ for (var i = 1; i <= ie; ++i)
6644N/A+ {
6644N/A+ tp <<= 4;
6644N/A+ if (tp <= ib && tp != BASE && i < ie)
6644N/A+ continue;
6644N/A+ z = z.multiply_integer (tp);
6644N/A+ tp = 1;
6644N/A+ }
6644N/A+ }
6644N/A+ }
6644N/A+
6644N/A+ re_sign = z.re_sign;
6644N/A+ im_sign = z.im_sign;
6644N/A+ re_exponent = z.re_exponent;
6644N/A+ im_exponent = z.im_exponent;
6644N/A+ for (var i = 0; i < z.re_fraction.length; i++)
6644N/A+ {
6644N/A+ re_fraction[i] = z.re_fraction[i];
6644N/A+ im_fraction[i] = z.im_fraction[i];
6644N/A+ }
6644N/A }
6644N/A
6644N/A public Number.complex (Number x, Number y)
6644N/A {
6644N/A- var tmp = MPFloat.init2 ((Precision) precision);
6644N/A- tmp.@set (x.re_num, Round.NEAREST);
6644N/A- re_num = tmp;
6644N/A- var tmp2 = MPFloat.init2 ((Precision) precision);
6644N/A- tmp2.@set (y.re_num, Round.NEAREST);
6644N/A- im_num = tmp2;
6644N/A+ re_sign = x.re_sign;
6644N/A+ re_exponent = x.re_exponent;
6644N/A+ for (var i = 0; i < im_fraction.length; i++)
6644N/A+ re_fraction[i] = x.re_fraction[i];
6644N/A+
6644N/A+ im_sign = y.re_sign;
6644N/A+ im_exponent = y.re_exponent;
6644N/A+ for (var i = 0; i < im_fraction.length; i++)
6644N/A+ im_fraction[i] = y.re_fraction[i];
6644N/A }
6644N/A
6644N/A public Number.polar (Number r, Number theta, AngleUnit unit = AngleUnit.RADIANS)
6644N/A@@ -125,35 +328,28 @@
6644N/A
6644N/A public Number.eulers ()
6644N/A {
6644N/A- var tmp = MPFloat.init2 ((Precision) precision);
6644N/A- var tmp2 = MPFloat.init2 ((Precision) precision);
6644N/A- tmp2.set_unsigned_integer (1, Round.NEAREST);
6644N/A- /* e^1, since mpfr doesn't have a function to return e */
6644N/A- tmp.exp (tmp2, Round.NEAREST);
6644N/A- re_num = tmp;
6644N/A- var tmp3 = MPFloat.init2 ((Precision) precision);
6644N/A- tmp3.set_unsigned_integer (0, Round.NEAREST);
6644N/A- im_num = tmp3;
6644N/A+ var z = new Number.integer (1).epowy ();
6644N/A+ re_sign = z.re_sign;
6644N/A+ im_sign = z.im_sign;
6644N/A+ re_exponent = z.re_exponent;
6644N/A+ im_exponent = z.im_exponent;
6644N/A+ for (var i = 0; i < z.re_fraction.length; i++)
6644N/A+ {
6644N/A+ re_fraction[i] = z.re_fraction[i];
6644N/A+ im_fraction[i] = z.im_fraction[i];
6644N/A+ }
6644N/A }
6644N/A
6644N/A public Number.i ()
6644N/A {
6644N/A- var tmp = MPFloat.init2 ((Precision) precision);
6644N/A- var tmp2 = MPFloat.init2 ((Precision) precision);
6644N/A- tmp.set_unsigned_integer (0, Round.NEAREST);
6644N/A- tmp2.set_unsigned_integer (1, Round.NEAREST);
6644N/A- re_num = tmp;
6644N/A- im_num = tmp2;
6644N/A+ im_sign = 1;
6644N/A+ im_exponent = 1;
6644N/A+ im_fraction[0] = 1;
6644N/A }
6644N/A
6644N/A public Number.pi ()
6644N/A {
6644N/A- var tmp = MPFloat.init2 ((Precision) precision);
6644N/A- tmp.const_pi (Round.NEAREST);
6644N/A- re_num = tmp;
6644N/A- var tmp2 = MPFloat.init2 ((Precision) precision);
6644N/A- tmp2.set_unsigned_integer (0, Round.NEAREST);
6644N/A- im_num = tmp2;
6644N/A+ Number.double (Math.PI);
6644N/A }
6644N/A
6644N/A /* Sets z to be a uniform random number in the range [0, 1] */
6644N/A@@ -162,42 +358,123 @@
6644N/A Number.double (Random.next_double ());
6644N/A }
6644N/A
6644N/A- ~Number ()
6644N/A- {
6644N/A- re_num.clear ();
6644N/A- im_num.clear ();
6644N/A- }
6644N/A-
6644N/A public int64 to_integer ()
6644N/A {
6644N/A- return re_num.get_signed_integer (Round.NEAREST);
6644N/A+ int64 z = 0;
6644N/A+
6644N/A+ /* |x| <= 1 */
6644N/A+ if (re_sign == 0 || re_exponent <= 0)
6644N/A+ return 0;
6644N/A+
6644N/A+ /* Multiply digits together */
6644N/A+ for (var i = 0; i < re_exponent; i++)
6644N/A+ {
6644N/A+ var t = z;
6644N/A+ z = z * BASE + re_fraction[i];
6644N/A+
6644N/A+ /* Check for overflow */
6644N/A+ if (z <= t)
6644N/A+ return 0;
6644N/A+ }
6644N/A+
6644N/A+ /* Validate result */
6644N/A+ var v = z;
6644N/A+ for (var i = re_exponent - 1; i >= 0; i--)
6644N/A+ {
6644N/A+ /* Get last digit */
6644N/A+ var digit = v - (v / BASE) * BASE;
6644N/A+ if (re_fraction[i] != digit)
6644N/A+ return 0;
6644N/A+
6644N/A+ v /= BASE;
6644N/A+ }
6644N/A+ if (v != 0)
6644N/A+ return 0;
6644N/A+
6644N/A+ return re_sign * z;
6644N/A }
6644N/A
6644N/A public uint64 to_unsigned_integer ()
6644N/A {
6644N/A- return re_num.get_unsigned_integer (Round.NEAREST);
6644N/A+ /* x <= 1 */
6644N/A+ if (re_sign <= 0 || re_exponent <= 0)
6644N/A+ return 0;
6644N/A+
6644N/A+ /* Multiply digits together */
6644N/A+ uint64 z = 0;
6644N/A+ for (var i = 0; i < re_exponent; i++)
6644N/A+ {
6644N/A+ var t = z;
6644N/A+ z = z * BASE + re_fraction[i];
6644N/A+
6644N/A+ /* Check for overflow */
6644N/A+ if (z <= t)
6644N/A+ return 0;
6644N/A+ }
6644N/A+
6644N/A+ /* Validate result */
6644N/A+ var v = z;
6644N/A+ for (var i = re_exponent - 1; i >= 0; i--)
6644N/A+ {
6644N/A+ /* Get last digit */
6644N/A+ var digit = (uint64) v - (v / BASE) * BASE;
6644N/A+ if (re_fraction[i] != digit)
6644N/A+ return 0;
6644N/A+
6644N/A+ v /= BASE;
6644N/A+ }
6644N/A+ if (v != 0)
6644N/A+ return 0;
6644N/A+
6644N/A+ return z;
6644N/A }
6644N/A
6644N/A public float to_float ()
6644N/A {
6644N/A- return re_num.get_float (Round.NEAREST);
6644N/A+ if (is_zero ())
6644N/A+ return 0f;
6644N/A+
6644N/A+ var z = 0f;
6644N/A+ for (var i = 0; i < T; i++)
6644N/A+ {
6644N/A+ if (re_fraction[i] != 0)
6644N/A+ z += re_fraction[i] * Math.powf (BASE, re_exponent - i - 1);
6644N/A+ }
6644N/A+
6644N/A+ if (re_sign < 0)
6644N/A+ return -z;
6644N/A+ else
6644N/A+ return z;
6644N/A }
6644N/A
6644N/A public double to_double ()
6644N/A {
6644N/A- return re_num.get_double (Round.NEAREST);
6644N/A+ if (is_zero ())
6644N/A+ return 0d;
6644N/A+
6644N/A+ var z = 0d;
6644N/A+ for (var i = 0; i < T; i++)
6644N/A+ {
6644N/A+ if (re_fraction[i] != 0)
6644N/A+ z += re_fraction[i] * Math.pow (BASE, re_exponent - i - 1);
6644N/A+ }
6644N/A+
6644N/A+ if (re_sign < 0)
6644N/A+ return -z;
6644N/A+ else
6644N/A+ return z;
6644N/A }
6644N/A
6644N/A /* Return true if the value is x == 0 */
6644N/A public bool is_zero ()
6644N/A {
6644N/A- return re_num.is_zero () && im_num.is_zero ();
6644N/A+ return re_sign == 0 && im_sign == 0;
6644N/A }
6644N/A
6644N/A /* Return true if x < 0 */
6644N/A public bool is_negative ()
6644N/A {
6644N/A- return re_num.sgn () < 0;
6644N/A+ return re_sign < 0;
6644N/A }
6644N/A
6644N/A /* Return true if x is integer */
6644N/A@@ -205,8 +482,35 @@
6644N/A {
6644N/A if (is_complex ())
6644N/A return false;
6644N/A+ /* Multiplication and division by 10000 is used to get around a
6644N/A+ * limitation to the "fix" for Sun bugtraq bug #4006391 in the
6644N/A+ * floor () routine in mp.c, when the re_exponent is less than 1.
6644N/A+ */
6644N/A+ var t3 = new Number.integer (10000);
6644N/A+ var t1 = multiply (t3);
6644N/A+ t1 = t1.divide (t3);
6644N/A+ var t2 = t1.floor ();
6644N/A+ return t1.equals (t2);
6644N/A
6644N/A- return re_num.is_integer () != 0;
6644N/A+ /* Correct way to check for integer */
6644N/A+ /*
6644N/A+
6644N/A+ // Zero is an integer
6644N/A+ if (is_zero ())
6644N/A+ return true;
6644N/A+
6644N/A+ // fractional
6644N/A+ if (re_exponent <= 0)
6644N/A+ return false;
6644N/A+
6644N/A+ // Look for fractional components
6644N/A+ for (var i = re_exponent; i < SIZE; i++)
6644N/A+ {
6644N/A+ if (re_fraction[i] != 0)
6644N/A+ return false;
6644N/A+ }
6644N/A+
6644N/A+ return true;*/
6644N/A }
6644N/A
6644N/A /* Return true if x is a positive integer */
6644N/A@@ -215,7 +519,7 @@
6644N/A if (is_complex ())
6644N/A return false;
6644N/A else
6644N/A- return re_num.sgn () >= 0 && is_integer ();
6644N/A+ return re_sign >= 0 && is_integer ();
6644N/A }
6644N/A
6644N/A /* Return true if x is a natural number (an integer ≥ 0) */
6644N/A@@ -224,34 +528,19 @@
6644N/A if (is_complex ())
6644N/A return false;
6644N/A else
6644N/A- return re_num.sgn () > 0 && is_integer ();
6644N/A+ return re_sign > 0 && is_integer ();
6644N/A }
6644N/A
6644N/A /* Return true if x has an imaginary component */
6644N/A public bool is_complex ()
6644N/A {
6644N/A- return !im_num.is_zero ();
6644N/A+ return im_sign != 0;
6644N/A }
6644N/A
6644N/A- /* Return error if overflow or underflow */
6644N/A- public static void check_flags ()
6644N/A- {
6644N/A- if (mpfr_is_underflow () != 0)
6644N/A- {
6644N/A- /* Translators: Error displayed when underflow error occured */
6644N/A- error = _("Underflow error");
6644N/A- }
6644N/A- else if (mpfr_is_overflow () != 0)
6644N/A- {
6644N/A- /* Translators: Error displayed when overflow error occured */
6644N/A- error = _("Overflow error");
6644N/A- }
6644N/A- }
6644N/A-
6644N/A /* Return true if x == y */
6644N/A public bool equals (Number y)
6644N/A {
6644N/A- return re_num.is_equal (y.re_num);
6644N/A+ return compare (y) == 0;
6644N/A }
6644N/A
6644N/A /* Returns:
6644N/A@@ -261,25 +550,62 @@
6644N/A */
6644N/A public int compare (Number y)
6644N/A {
6644N/A- return re_num.cmp (y.re_num);
6644N/A+ if (re_sign != y.re_sign)
6644N/A+ {
6644N/A+ if (re_sign > y.re_sign)
6644N/A+ return 1;
6644N/A+ else
6644N/A+ return -1;
6644N/A+ }
6644N/A+
6644N/A+ /* x = y = 0 */
6644N/A+ if (is_zero ())
6644N/A+ return 0;
6644N/A+
6644N/A+ /* See if numbers are of different magnitude */
6644N/A+ if (re_exponent != y.re_exponent)
6644N/A+ {
6644N/A+ if (re_exponent > y.re_exponent)
6644N/A+ return re_sign;
6644N/A+ else
6644N/A+ return -re_sign;
6644N/A+ }
6644N/A+
6644N/A+ /* Compare fractions */
6644N/A+ for (var i = 0; i < SIZE; i++)
6644N/A+ {
6644N/A+ if (re_fraction[i] == y.re_fraction[i])
6644N/A+ continue;
6644N/A+
6644N/A+ if (re_fraction[i] > y.re_fraction[i])
6644N/A+ return re_sign;
6644N/A+ else
6644N/A+ return -re_sign;
6644N/A+ }
6644N/A+
6644N/A+ /* x = y */
6644N/A+ return 0;
6644N/A }
6644N/A
6644N/A /* Sets z = sgn (x) */
6644N/A public Number sgn ()
6644N/A {
6644N/A- var z = new Number.integer (re_num.sgn ());
6644N/A- return z;
6644N/A+ if (is_zero ())
6644N/A+ return new Number.integer (0);
6644N/A+ else if (is_negative ())
6644N/A+ return new Number.integer (-1);
6644N/A+ else
6644N/A+ return new Number.integer (1);
6644N/A }
6644N/A
6644N/A /* Sets z = −x */
6644N/A public Number invert_sign ()
6644N/A {
6644N/A- var tmp = MPFloat.init2 ((Precision) precision);
6644N/A- tmp.neg (re_num, Round.NEAREST);
6644N/A- var z = new Number.mpfloat (tmp);
6644N/A- var tmp_im = z.im_num;
6644N/A- tmp_im.neg (im_num, Round.NEAREST);
6644N/A- z.im_num = tmp_im;
6644N/A+ var z = copy ();
6644N/A+
6644N/A+ z.re_sign = -z.re_sign;
6644N/A+ z.im_sign = -z.im_sign;
6644N/A+
6644N/A return z;
6644N/A }
6644N/A
6644N/A@@ -294,13 +620,13 @@
6644N/A x_real = x_real.multiply (x_real);
6644N/A x_im = x_im.multiply (x_im);
6644N/A var z = x_real.add (x_im);
6644N/A- return z.sqrt ();
6644N/A+ return z.root (2);
6644N/A }
6644N/A else
6644N/A {
6644N/A- var tmp = MPFloat.init2 ((Precision) precision);
6644N/A- tmp.abs (re_num, Round.NEAREST);
6644N/A- var z = new Number.mpfloat (tmp);
6644N/A+ var z = copy ();
6644N/A+ if (z.re_sign < 0)
6644N/A+ z.re_sign = -z.re_sign;
6644N/A return z;
6644N/A }
6644N/A }
6644N/A@@ -311,7 +637,7 @@
6644N/A if (is_zero ())
6644N/A {
6644N/A /* Translators: Error display when attempting to take argument of zero */
6644N/A- error = _("Argument not defined for zero");
6644N/A+ mperr (_("Argument not defined for zero"));
6644N/A return new Number.integer (0);
6644N/A }
6644N/A
6644N/A@@ -355,12 +681,8 @@
6644N/A /* Sets z = ‾̅x */
6644N/A public Number conjugate ()
6644N/A {
6644N/A- var tmp = MPFloat.init2 ((Precision) precision);
6644N/A- tmp.neg (im_num, Round.NEAREST);
6644N/A var z = copy ();
6644N/A- var tmp2 = z.im_num;
6644N/A- tmp2.clear ();
6644N/A- z.im_num = tmp;
6644N/A+ z.im_sign = -z.im_sign;
6644N/A return z;
6644N/A }
6644N/A
6644N/A@@ -368,11 +690,13 @@
6644N/A public Number real_component ()
6644N/A {
6644N/A var z = copy ();
6644N/A- var tmp = z.im_num;
6644N/A- tmp.clear ();
6644N/A- tmp = MPFloat.init2 ((Precision) precision);
6644N/A- tmp.set_unsigned_integer (0, Round.NEAREST);
6644N/A- z.im_num = tmp;
6644N/A+
6644N/A+ /* Clear imaginary component */
6644N/A+ z.im_sign = 0;
6644N/A+ z.im_exponent = 0;
6644N/A+ for (var i = 0; i < z.im_fraction.length; i++)
6644N/A+ z.im_fraction[i] = 0;
6644N/A+
6644N/A return z;
6644N/A }
6644N/A
6644N/A@@ -380,30 +704,65 @@
6644N/A public Number imaginary_component ()
6644N/A {
6644N/A /* Copy imaginary component to real component */
6644N/A- var z = copy ();
6644N/A- var tmp = z.re_num;
6644N/A- tmp.clear ();
6644N/A- z.re_num = z.im_num;
6644N/A- tmp = MPFloat.init2 ((Precision) precision);
6644N/A- tmp.set_unsigned_integer (0, Round.NEAREST);
6644N/A- z.im_num = tmp;
6644N/A+ var z = new Number ();
6644N/A+ z.re_sign = im_sign;
6644N/A+ z.re_exponent = im_exponent;
6644N/A+
6644N/A+ for (var i = 0; i < z.im_fraction.length; i++)
6644N/A+ z.re_fraction[i] = im_fraction[i];
6644N/A+
6644N/A+ /* Clear (old) imaginary component */
6644N/A+ z.im_sign = 0;
6644N/A+ z.im_exponent = 0;
6644N/A+ for (var i = 0; i < z.im_fraction.length; i++)
6644N/A+ z.im_fraction[i] = 0;
6644N/A+
6644N/A return z;
6644N/A }
6644N/A
6644N/A public Number integer_component ()
6644N/A {
6644N/A- var tmp = MPFloat.init2 ((Precision) precision);
6644N/A- tmp.trunc (re_num);
6644N/A- var z = new Number.mpfloat (tmp);
6644N/A+ /* Clear re_fraction */
6644N/A+ var z = copy ();
6644N/A+ for (var i = z.re_exponent; i < SIZE; i++)
6644N/A+ z.re_fraction[i] = 0;
6644N/A+ z.im_sign = 0;
6644N/A+ z.im_exponent = 0;
6644N/A+ for (var i = 0; i < z.im_fraction.length; i++)
6644N/A+ z.im_fraction[i] = 0;
6644N/A+
6644N/A return z;
6644N/A+
6644N/A }
6644N/A
6644N/A /* Sets z = x mod 1 */
6644N/A public Number fractional_component ()
6644N/A {
6644N/A- var tmp = MPFloat.init2 ((Precision) precision);
6644N/A- tmp.frac (re_num, Round.NEAREST);
6644N/A- var z = new Number.mpfloat (tmp);
6644N/A+ /* fractional component of zero is 0 */
6644N/A+ if (is_zero ())
6644N/A+ return new Number.integer (0);
6644N/A+
6644N/A+ /* All fractional */
6644N/A+ if (re_exponent <= 0)
6644N/A+ return this;
6644N/A+
6644N/A+ /* Shift fractional component */
6644N/A+ var shift = re_exponent;
6644N/A+ for (var i = shift; i < SIZE && re_fraction[i] == 0; i++)
6644N/A+ shift++;
6644N/A+ var z = new Number.integer (0);
6644N/A+ z.re_sign = re_sign;
6644N/A+ z.re_exponent = re_exponent - shift;
6644N/A+ for (var i = 0; i < SIZE; i++)
6644N/A+ {
6644N/A+ if (i + shift >= SIZE)
6644N/A+ z.re_fraction[i] = 0;
6644N/A+ else
6644N/A+ z.re_fraction[i] = re_fraction[i + shift];
6644N/A+ }
6644N/A+ if (z.re_fraction[0] == 0)
6644N/A+ z.re_sign = 0;
6644N/A+
6644N/A return z;
6644N/A }
6644N/A
6644N/A@@ -416,9 +775,36 @@
6644N/A /* Sets z = ⌊x⌋ */
6644N/A public Number floor ()
6644N/A {
6644N/A- var tmp = MPFloat.init2 ((Precision) precision);
6644N/A- tmp.floor (re_num);
6644N/A- var z = new Number.mpfloat (tmp);
6644N/A+ /* Integer component of zero = 0 */
6644N/A+ if (is_zero ())
6644N/A+ return this;
6644N/A+
6644N/A+ /* If all fractional then no integer component */
6644N/A+ if (re_exponent <= 0)
6644N/A+ {
6644N/A+ if (is_negative ())
6644N/A+ return new Number.integer (-1);
6644N/A+ else
6644N/A+ return new Number.integer (0);
6644N/A+ }
6644N/A+
6644N/A+ /* Clear fractional component */
6644N/A+ var z = copy ();
6644N/A+ var have_fraction = false;
6644N/A+ for (var i = z.re_exponent; i < SIZE; i++)
6644N/A+ {
6644N/A+ if (z.re_fraction[i] != 0)
6644N/A+ have_fraction = true;
6644N/A+ z.re_fraction[i] = 0;
6644N/A+ }
6644N/A+ z.im_sign = 0;
6644N/A+ z.im_exponent = 0;
6644N/A+ for (var i = 0; i < z.im_fraction.length; i++)
6644N/A+ z.im_fraction[i] = 0;
6644N/A+
6644N/A+ if (have_fraction && is_negative ())
6644N/A+ z = z.add (new Number.integer (-1));
6644N/A+
6644N/A return z;
6644N/A }
6644N/A
6644N/A@@ -425,19 +811,29 @@
6644N/A /* Sets z = ⌈x⌉ */
6644N/A public Number ceiling ()
6644N/A {
6644N/A- var tmp = MPFloat.init2 ((Precision) precision);
6644N/A- tmp.ceil (re_num);
6644N/A- var z = new Number.mpfloat (tmp);
6644N/A- return z;
6644N/A+ var z = floor ();
6644N/A+ var f = fractional_component ();
6644N/A+ if (f.is_zero ())
6644N/A+ return z;
6644N/A+ return z.add (new Number.integer (1));
6644N/A }
6644N/A
6644N/A /* Sets z = [x] */
6644N/A public Number round ()
6644N/A {
6644N/A- var tmp = MPFloat.init2 ((Precision) precision);
6644N/A- tmp.round (re_num);
6644N/A- var z = new Number.mpfloat (tmp);
6644N/A- return z;
6644N/A+ var do_floor = !is_negative ();
6644N/A+
6644N/A+ var f = fractional_component ();
6644N/A+ f = f.multiply_integer (2);
6644N/A+ f = f.abs ();
6644N/A+ if (f.compare (new Number.integer (1)) >= 0)
6644N/A+ do_floor = !do_floor;
6644N/A+
6644N/A+ if (do_floor)
6644N/A+ return floor ();
6644N/A+ else
6644N/A+ return ceiling ();
6644N/A+
6644N/A }
6644N/A
6644N/A /* Sets z = 1 ÷ x */
6644N/A@@ -482,24 +878,10 @@
6644N/A /* Sets z = x^y */
6644N/A public Number xpowy (Number y)
6644N/A {
6644N/A- /* 0^-n invalid */
6644N/A- if (is_zero () && y.is_negative ())
6644N/A+ if (y.is_integer ())
6644N/A+ return xpowy_integer (y.to_integer ());
6644N/A+ else
6644N/A {
6644N/A- /* Translators: Error displayed when attempted to raise 0 to a negative re_exponent */
6644N/A- error = _("The power of zero is undefined for a negative exponent");
6644N/A- return new Number.integer (0);
6644N/A- }
6644N/A-
6644N/A- /* 0^0 is indeterminate */
6644N/A- if (is_zero () && y.is_zero ())
6644N/A- {
6644N/A- /* Translators: Error displayed when attempted to raise 0 to power of zero */
6644N/A- error = _("Zero raised to zero is undefined");
6644N/A- return new Number.integer (0);
6644N/A- }
6644N/A-
6644N/A- if (!y.is_integer ())
6644N/A- {
6644N/A var reciprocal = y.reciprocal ();
6644N/A if (reciprocal.is_integer ())
6644N/A return root (reciprocal.to_integer ());
6644N/A@@ -506,29 +888,6 @@
6644N/A else
6644N/A return pwr (y);
6644N/A }
6644N/A-
6644N/A- Number t;
6644N/A- Number t2;
6644N/A- if (y.is_negative ())
6644N/A- {
6644N/A- t = reciprocal ();
6644N/A- t2 = y.invert_sign ();
6644N/A- }
6644N/A- else
6644N/A- {
6644N/A- t = this;
6644N/A- t2 = y;
6644N/A- }
6644N/A-
6644N/A- var tmp = MPFloat.init2 ((Precision) precision);
6644N/A- tmp.power (t.re_num, t2.re_num, Round.NEAREST);
6644N/A- var z = new Number.mpfloat (tmp);
6644N/A- var tmp2 = z.im_num;
6644N/A- tmp2.clear ();
6644N/A- tmp = MPFloat.init2 ((Precision) precision);
6644N/A- tmp.@set (t.im_num, Round.NEAREST);
6644N/A- z.im_num = tmp;
6644N/A- return z;
6644N/A }
6644N/A
6644N/A /* Sets z = x^y */
6644N/A@@ -538,7 +897,7 @@
6644N/A if (is_zero () && n < 0)
6644N/A {
6644N/A /* Translators: Error displayed when attempted to raise 0 to a negative re_exponent */
6644N/A- error = _("The power of zero is undefined for a negative exponent");
6644N/A+ mperr (_("The power of zero is undefined for a negative exponent"));
6644N/A return new Number.integer (0);
6644N/A }
6644N/A
6644N/A@@ -546,10 +905,22 @@
6644N/A if (is_zero () && n == 0)
6644N/A {
6644N/A /* Translators: Error displayed when attempted to raise 0 to power of zero */
6644N/A- error = _("Zero raised to zero is undefined");
6644N/A+ mperr (_("Zero raised to zero is undefined"));
6644N/A return new Number.integer (0);
6644N/A }
6644N/A
6644N/A+ /* x^0 = 1 */
6644N/A+ if (n == 0)
6644N/A+ return new Number.integer (1);
6644N/A+
6644N/A+ /* 0^n = 0 */
6644N/A+ if (is_zero ())
6644N/A+ return new Number.integer (0);
6644N/A+
6644N/A+ /* x^1 = x */
6644N/A+ if (n == 1)
6644N/A+ return this;
6644N/A+
6644N/A Number t;
6644N/A if (n < 0)
6644N/A {
6644N/A@@ -557,44 +928,17 @@
6644N/A n = -n;
6644N/A }
6644N/A else
6644N/A- {
6644N/A t = this;
6644N/A- }
6644N/A
6644N/A- var tmp = MPFloat.init2 ((Precision) precision);
6644N/A- tmp.power_signed_integer (t.re_num, (long) n, Round.NEAREST);
6644N/A- var z = new Number.mpfloat (tmp);
6644N/A- var tmp2 = z.im_num;
6644N/A- tmp2.clear ();
6644N/A- tmp = MPFloat.init2 ((Precision) precision);
6644N/A- tmp.@set (t.im_num, Round.NEAREST);
6644N/A- z.im_num = tmp;
6644N/A- return z;
6644N/A- }
6644N/A-
6644N/A- private Number pwr (Number y)
6644N/A- {
6644N/A- /* (-x)^y imaginary */
6644N/A- /* FIXME: Make complex numbers optional */
6644N/A- /*if (re_sign < 0)
6644N/A+ var z = new Number.integer (1);
6644N/A+ while (n != 0)
6644N/A {
6644N/A- mperr (_("The power of negative numbers is only defined for integer exponents"));
6644N/A- return new Number.integer (0);
6644N/A- }*/
6644N/A-
6644N/A- /* 0^y = 0, 0^-y undefined */
6644N/A- if (is_zero ())
6644N/A- {
6644N/A- if (y.is_negative ())
6644N/A- error = _("The power of zero is undefined for a negative exponent");
6644N/A- return new Number.integer (0);
6644N/A+ if (n % 2 == 1)
6644N/A+ z = z.multiply (t);
6644N/A+ t = t.multiply (t);
6644N/A+ n = n / 2;
6644N/A }
6644N/A-
6644N/A- /* x^0 = 1 */
6644N/A- if (y.is_zero ())
6644N/A- return new Number.integer (1);
6644N/A-
6644N/A- return y.multiply (ln ()).epowy ();
6644N/A+ return z;
6644N/A }
6644N/A
6644N/A /* Sets z = n√x */
6644N/A@@ -623,10 +967,22 @@
6644N/A /* Sets z = √x */
6644N/A public Number sqrt ()
6644N/A {
6644N/A- var tmp = MPFloat.init2 ((Precision) precision);
6644N/A- tmp.sqrt (re_num, Round.NEAREST);
6644N/A- var z = new Number.mpfloat (tmp);
6644N/A- return z;
6644N/A+ if (is_zero ())
6644N/A+ return this;
6644N/A+
6644N/A+ /* FIXME: Make complex numbers optional */
6644N/A+ /*if (re_sign < 0)
6644N/A+ {
6644N/A+ mperr (_("Square root is undefined for negative values"));
6644N/A+ return new Number.integer (0);
6644N/A+ }*/
6644N/A+ else
6644N/A+ {
6644N/A+ var t = root (-2);
6644N/A+ var z = multiply (t);
6644N/A+ return z.ext (t.re_fraction[0], z.re_fraction[0]);
6644N/A+ }
6644N/A+
6644N/A }
6644N/A
6644N/A /* Sets z = ln x */
6644N/A@@ -636,7 +992,7 @@
6644N/A if (is_zero ())
6644N/A {
6644N/A /* Translators: Error displayed when attempting to take logarithm of zero */
6644N/A- error = _("Logarithm of zero is undefined");
6644N/A+ mperr (_("Logarithm of zero is undefined"));
6644N/A return new Number.integer (0);
6644N/A }
6644N/A
6644N/A@@ -669,7 +1025,7 @@
6644N/A if (is_zero ())
6644N/A {
6644N/A /* Translators: Error displayed when attempting to take logarithm of zero */
6644N/A- error = _("Logarithm of zero is undefined");
6644N/A+ mperr (_("Logarithm of zero is undefined"));
6644N/A return new Number.integer (0);
6644N/A }
6644N/A
6644N/A@@ -691,17 +1047,17 @@
6644N/A if (is_negative () || is_complex ())
6644N/A {
6644N/A /* Translators: Error displayed when attempted take the factorial of a negative or complex number */
6644N/A- error = _("Factorial is only defined for non-negative real numbers");
6644N/A+ mperr (_("Factorial is only defined for non-negative real numbers"));
6644N/A return new Number.integer (0);
6644N/A }
6644N/A
6644N/A- var tmp = add (new Number.integer (1));
6644N/A- var tmp2 = MPFloat.init2 ((Precision) precision);
6644N/A+ var val = to_double ();
6644N/A
6644N/A /* Factorial(x) = Gamma(x+1) - This is the formula used to calculate Factorial.*/
6644N/A- tmp2.gamma (tmp.re_num, Round.NEAREST);
6644N/A+ var fact = Math.tgamma(val+1);
6644N/A
6644N/A- return new Number.mpfloat (tmp2);
6644N/A+ return new Number.double(fact);
6644N/A+
6644N/A }
6644N/A
6644N/A /* Convert to integer - if couldn't be converted then the factorial would be too big anyway */
6644N/A@@ -716,41 +1072,23 @@
6644N/A /* Sets z = x + y */
6644N/A public Number add (Number y)
6644N/A {
6644N/A- if (is_complex () || y.is_complex ())
6644N/A- {
6644N/A- Number real_z, im_z;
6644N/A-
6644N/A- var real_x = real_component ();
6644N/A- var im_x = imaginary_component ();
6644N/A- var real_y = y.real_component ();
6644N/A- var im_y = y.imaginary_component ();
6644N/A-
6644N/A- real_z = real_x.add_real (real_y);
6644N/A- im_z = im_x.add_real (im_y);
6644N/A-
6644N/A- return new Number.complex (real_z, im_z);
6644N/A- }
6644N/A- else
6644N/A- return add_real (y);
6644N/A+ return add_with_sign (1, y);
6644N/A }
6644N/A
6644N/A- public Number add_real (Number y)
6644N/A- {
6644N/A- var tmp = MPFloat.init2 ((Precision) precision);
6644N/A- tmp.add (re_num, y.re_num, Round.NEAREST);
6644N/A- var z = new Number.mpfloat (tmp);
6644N/A- return z;
6644N/A- }
6644N/A-
6644N/A /* Sets z = x − y */
6644N/A public Number subtract (Number y)
6644N/A {
6644N/A- return add (y.invert_sign ());
6644N/A+ return add_with_sign (-1, y);
6644N/A }
6644N/A
6644N/A /* Sets z = x × y */
6644N/A public Number multiply (Number y)
6644N/A {
6644N/A+ /* x*0 = 0*y = 0 */
6644N/A+ if (is_zero () || y.is_zero ())
6644N/A+ return new Number.integer (0);
6644N/A+
6644N/A+ /* (a+bi)(c+di) = (ac-bd)+(ad+bc)i */
6644N/A if (is_complex () || y.is_complex ())
6644N/A {
6644N/A Number t1, t2, real_z, im_z;
6644N/A@@ -774,23 +1112,17 @@
6644N/A return multiply_real (y);
6644N/A }
6644N/A
6644N/A- public Number multiply_real (Number y)
6644N/A- {
6644N/A- var z = new Number.integer (0);
6644N/A- var tmp = z.re_num;
6644N/A- tmp.multiply (re_num, y.re_num, Round.NEAREST);
6644N/A- z.re_num = tmp;
6644N/A- return z;
6644N/A- }
6644N/A-
6644N/A /* Sets z = x × y */
6644N/A public Number multiply_integer (int64 y)
6644N/A {
6644N/A- var z = new Number.integer (0);
6644N/A- var tmp = z.re_num;
6644N/A- tmp.multiply_signed_integer (re_num, (long) y, Round.NEAREST);
6644N/A- z.re_num = tmp;
6644N/A- return z;
6644N/A+ if (is_complex ())
6644N/A+ {
6644N/A+ var re_z = real_component ().multiply_integer_real (y);;
6644N/A+ var im_z = imaginary_component ().multiply_integer_real (y);
6644N/A+ return new Number.complex (re_z, im_z);
6644N/A+ }
6644N/A+ else
6644N/A+ return multiply_integer_real (y);
6644N/A }
6644N/A
6644N/A /* Sets z = x ÷ y */
6644N/A@@ -799,31 +1131,23 @@
6644N/A if (y.is_zero ())
6644N/A {
6644N/A /* Translators: Error displayed attempted to divide by zero */
6644N/A- error = _("Division by zero is undefined");
6644N/A+ mperr (_("Division by zero is undefined"));
6644N/A return new Number.integer (0);
6644N/A }
6644N/A+ /* 0/y = 0 */
6644N/A+ if (is_zero ())
6644N/A+ return this;
6644N/A
6644N/A- if (is_complex () || y.is_complex ())
6644N/A- {
6644N/A- var a = real_component ();
6644N/A- var b = imaginary_component ();
6644N/A- var c = y.real_component ();
6644N/A- var d = y.imaginary_component ();
6644N/A+ /* z = x × y⁻¹ */
6644N/A+ /* FIXME: Set re_exponent to zero to avoid overflow in multiply??? */
6644N/A+ var t = y.reciprocal ();
6644N/A+ var ie = t.re_exponent;
6644N/A+ t.re_exponent = 0;
6644N/A+ var i = t.re_fraction[0];
6644N/A+ var z = multiply (t);
6644N/A+ z = z.ext (i, z.re_fraction[0]);
6644N/A+ z.re_exponent += ie;
6644N/A
6644N/A- var tmp = a.multiply (c).add (b.multiply (d));
6644N/A- var tmp_2 = c.xpowy_integer (2).add (d.xpowy_integer (2));
6644N/A- var z_1 = tmp.divide (tmp_2);
6644N/A-
6644N/A- tmp = b.multiply (c).subtract (a.multiply (d));
6644N/A- var z_2 = tmp.divide (tmp_2);
6644N/A-
6644N/A- var z = new Number.complex (z_1, z_2);
6644N/A- return z;
6644N/A- }
6644N/A-
6644N/A- var tmp = MPFloat.init2 ((Precision) precision);
6644N/A- tmp.divide (re_num, y.re_num, Round.NEAREST);
6644N/A- var z = new Number.mpfloat (tmp);
6644N/A return z;
6644N/A }
6644N/A
6644N/A@@ -830,7 +1154,14 @@
6644N/A /* Sets z = x ÷ y */
6644N/A public Number divide_integer (int64 y)
6644N/A {
6644N/A- return divide (new Number.integer (y));
6644N/A+ if (is_complex ())
6644N/A+ {
6644N/A+ var re_z = real_component ().divide_integer_real (y);
6644N/A+ var im_z = imaginary_component ().divide_integer_real (y);
6644N/A+ return new Number.complex (re_z, im_z);
6644N/A+ }
6644N/A+ else
6644N/A+ return divide_integer_real (y);
6644N/A }
6644N/A
6644N/A /* Sets z = x mod y */
6644N/A@@ -839,7 +1170,7 @@
6644N/A if (!is_integer () || !y.is_integer ())
6644N/A {
6644N/A /* Translators: Error displayed when attemping to do a modulus division on non-integer numbers */
6644N/A- error = _("Modulus division is only defined for integers");
6644N/A+ mperr (_("Modulus division is only defined for integers"));
6644N/A return new Number.integer (0);
6644N/A }
6644N/A
6644N/A@@ -857,7 +1188,7 @@
6644N/A /* Sets z = x ^ y mod p */
6644N/A public Number modular_exponentiation (Number exp, Number mod)
6644N/A {
6644N/A- var base_value = copy ();
6644N/A+ var base_value = this.copy ();
6644N/A if (exp.is_negative ())
6644N/A base_value = base_value.reciprocal ();
6644N/A var exp_value = exp.abs ();
6644N/A@@ -927,51 +1258,87 @@
6644N/A public Number tan (AngleUnit unit = AngleUnit.RADIANS)
6644N/A {
6644N/A /* Check for undefined values */
6644N/A- var x_radians = to_radians (unit);
6644N/A- var check = x_radians.subtract (new Number.pi ().divide_integer (2)).divide (new Number.pi ());
6644N/A-
6644N/A- if (check.is_integer ())
6644N/A+ var cos_x = cos (unit);
6644N/A+ if (cos_x.is_zero ())
6644N/A {
6644N/A /* Translators: Error displayed when tangent value is undefined */
6644N/A- error = _("Tangent is undefined for angles that are multiples of π (180°) from π∕2 (90°)");
6644N/A+ mperr (_("Tangent is undefined for angles that are multiples of π (180°) from π∕2 (90°)"));
6644N/A return new Number.integer (0);
6644N/A }
6644N/A
6644N/A- var tmp = MPFloat.init2 ((Precision) precision);
6644N/A- tmp.tan (x_radians.re_num, Round.NEAREST);
6644N/A- var z = new Number.mpfloat (tmp);
6644N/A- return z;
6644N/A+ /* tan (x) = sin (x) / cos (x) */
6644N/A+ return sin (unit).divide (cos_x);
6644N/A }
6644N/A
6644N/A /* Sets z = sin⁻¹ x */
6644N/A public Number asin (AngleUnit unit = AngleUnit.RADIANS)
6644N/A {
6644N/A- if (compare (new Number.integer (1)) > 0 || compare (new Number.integer (-1)) < 0)
6644N/A- {
6644N/A- /* Translators: Error displayed when inverse sine value is undefined */
6644N/A- error = _("Inverse sine is undefined for values outside [-1, 1]");
6644N/A+ /* asin⁻¹(0) = 0 */
6644N/A+ if (is_zero ())
6644N/A return new Number.integer (0);
6644N/A+
6644N/A+ /* sin⁻¹(x) = tan⁻¹(x / √(1 - x²)), |x| < 1 */
6644N/A+ if (re_exponent <= 0)
6644N/A+ {
6644N/A+ var t1 = new Number.integer (1);
6644N/A+ var t2 = t1;
6644N/A+ t1 = t1.subtract (this);
6644N/A+ t2 = t2.add (this);
6644N/A+ t2 = t1.multiply (t2);
6644N/A+ t2 = t2.root (-2);
6644N/A+ var z = multiply (t2);
6644N/A+ z = z.atan (unit);
6644N/A+ return z;
6644N/A }
6644N/A
6644N/A- var tmp = MPFloat.init2 ((Precision) precision);
6644N/A- tmp.asin (re_num, Round.NEAREST);
6644N/A- var z = new Number.mpfloat (tmp);
6644N/A- return z.from_radians (unit);
6644N/A+ /* sin⁻¹(1) = π/2, sin⁻¹(-1) = -π/2 */
6644N/A+ var t2 = new Number.integer (re_sign);
6644N/A+ if (equals (t2))
6644N/A+ {
6644N/A+ var z = new Number.pi ().divide_integer (2 * t2.re_sign);
6644N/A+ return z.from_radians (unit);
6644N/A+ }
6644N/A+
6644N/A+ /* Translators: Error displayed when inverse sine value is undefined */
6644N/A+ mperr (_("Inverse sine is undefined for values outside [-1, 1]"));
6644N/A+ return new Number.integer (0);
6644N/A }
6644N/A
6644N/A /* Sets z = cos⁻¹ x */
6644N/A public Number acos (AngleUnit unit = AngleUnit.RADIANS)
6644N/A {
6644N/A- if (compare (new Number.integer (1)) > 0 || compare (new Number.integer (-1)) < 0)
6644N/A+ var pi = new Number.pi ();
6644N/A+ var t1 = new Number.integer (1);
6644N/A+ var n1 = new Number.integer (-1);
6644N/A+
6644N/A+ Number z;
6644N/A+ if (compare (t1) > 0 || compare (n1) < 0)
6644N/A {
6644N/A /* Translators: Error displayed when inverse cosine value is undefined */
6644N/A- error = _("Inverse cosine is undefined for values outside [-1, 1]");
6644N/A- return new Number.integer (0);
6644N/A+ mperr (_("Inverse cosine is undefined for values outside [-1, 1]"));
6644N/A+ z = new Number.integer (0);
6644N/A }
6644N/A+ else if (is_zero ())
6644N/A+ z = pi.divide_integer (2);
6644N/A+ else if (equals (t1))
6644N/A+ z = new Number.integer (0);
6644N/A+ else if (equals (n1))
6644N/A+ z = pi;
6644N/A+ else
6644N/A+ {
6644N/A+ /* cos⁻¹(x) = tan⁻¹(√(1 - x²) / x) */
6644N/A+ Number y;
6644N/A+ var t2 = multiply (this);
6644N/A+ t2 = t1.subtract (t2);
6644N/A+ t2 = t2.sqrt ();
6644N/A+ t2 = t2.divide (this);
6644N/A+ y = t2.atan (AngleUnit.RADIANS);
6644N/A+ if (re_sign > 0)
6644N/A+ z = y;
6644N/A+ else
6644N/A+ z = y.add (pi);
6644N/A+ }
6644N/A
6644N/A- var tmp = MPFloat.init2 ((Precision) precision);
6644N/A- tmp.acos (re_num, Round.NEAREST);
6644N/A- var z = new Number.mpfloat (tmp);
6644N/A return z.from_radians (unit);
6644N/A }
6644N/A
6644N/A@@ -978,9 +1345,63 @@
6644N/A /* Sets z = tan⁻¹ x */
6644N/A public Number atan (AngleUnit unit = AngleUnit.RADIANS)
6644N/A {
6644N/A- var tmp = MPFloat.init2 ((Precision) precision);
6644N/A- tmp.atan (re_num, Round.NEAREST);
6644N/A- var z = new Number.mpfloat (tmp);
6644N/A+ if (is_zero ())
6644N/A+ return new Number.integer (0);
6644N/A+
6644N/A+ var t2 = this;
6644N/A+ var rx = 0f;
6644N/A+ if (re_exponent.abs () <= 2)
6644N/A+ rx = to_float ();
6644N/A+
6644N/A+ /* REDUCE ARGUMENT IF NECESSARY BEFORE USING SERIES */
6644N/A+ var q = 1;
6644N/A+ Number z;
6644N/A+ while (t2.re_exponent >= 0)
6644N/A+ {
6644N/A+ if (t2.re_exponent == 0 && 2 * (t2.re_fraction[0] + 1) <= BASE)
6644N/A+ break;
6644N/A+
6644N/A+ q *= 2;
6644N/A+
6644N/A+ /* t = t / (√(t² + 1) + 1) */
6644N/A+ z = t2.multiply (t2);
6644N/A+ z = z.add (new Number.integer (1));
6644N/A+ z = z.sqrt ();
6644N/A+ z = z.add (new Number.integer (1));
6644N/A+ t2 = t2.divide (z);
6644N/A+ }
6644N/A+
6644N/A+ /* USE POWER SERIES NOW ARGUMENT IN (-0.5, 0.5) */
6644N/A+ z = t2;
6644N/A+ var t1 = t2.multiply (t2);
6644N/A+
6644N/A+ /* SERIES LOOP. REDUCE T IF POSSIBLE. */
6644N/A+ for (var i = 1; ; i += 2)
6644N/A+ {
6644N/A+ if (T + 2 + t2.re_exponent <= 1)
6644N/A+ break;
6644N/A+
6644N/A+ t2 = t2.multiply (t1).multiply_integer (-i).divide_integer (i + 2);
6644N/A+
6644N/A+ z = z.add (t2);
6644N/A+ if (t2.is_zero ())
6644N/A+ break;
6644N/A+ }
6644N/A+
6644N/A+ /* CORRECT FOR ARGUMENT REDUCTION */
6644N/A+ z = z.multiply_integer (q);
6644N/A+
6644N/A+ /* CHECK THAT RELATIVE ERROR LESS THAN 0.01 UNLESS re_exponent
6644N/A+ * OF X IS LARGE (WHEN ATAN MIGHT NOT WORK)
6644N/A+ */
6644N/A+ if (re_exponent.abs () <= 2)
6644N/A+ {
6644N/A+ float ry = z.to_float ();
6644N/A+ /* THE FOLLOWING MESSAGE MAY INDICATE THAT B**(T-1) IS TOO SMALL. */
6644N/A+ if (Math.fabs (ry - Math.atan (rx)) >= Math.fabs (ry) * 0.01)
6644N/A+ mperr ("*** ERROR OCCURRED IN ATAN, RESULT INCORRECT ***");
6644N/A+ }
6644N/A+
6644N/A return z.from_radians (unit);
6644N/A }
6644N/A
6644N/A@@ -987,27 +1408,88 @@
6644N/A /* Sets z = sinh x */
6644N/A public Number sinh ()
6644N/A {
6644N/A- var tmp = MPFloat.init2 ((Precision) precision);
6644N/A- tmp.sinh (re_num, Round.NEAREST);
6644N/A- var z = new Number.mpfloat (tmp);
6644N/A- return z;
6644N/A+ /* sinh (0) = 0 */
6644N/A+ if (is_zero ())
6644N/A+ return new Number.integer (0);
6644N/A+
6644N/A+ /* WORK WITH ABS (X) */
6644N/A+ var abs_x = abs ();
6644N/A+
6644N/A+ /* If |x| < 1 USE EXP TO AVOID CANCELLATION, otherwise IF TOO LARGE EPOWY GIVES ERROR MESSAGE */
6644N/A+ Number z;
6644N/A+ if (abs_x.re_exponent <= 0)
6644N/A+ {
6644N/A+ /* ((e^|x| + 1) * (e^|x| - 1)) / e^|x| */
6644N/A+ // FIXME: Solves to e^|x| - e^-|x|, why not lower branch always? */
6644N/A+ var exp_x = abs_x.epowy ();
6644N/A+ var a = exp_x.add (new Number.integer (1));
6644N/A+ var b = exp_x.add (new Number.integer (-1));
6644N/A+ z = a.multiply (b);
6644N/A+ z = z.divide (exp_x);
6644N/A+ }
6644N/A+ else
6644N/A+ {
6644N/A+ /* e^|x| - e^-|x| */
6644N/A+ var exp_x = abs_x.epowy ();
6644N/A+ z = exp_x.reciprocal ();
6644N/A+ z = exp_x.subtract (z);
6644N/A+ }
6644N/A+
6644N/A+ /* DIVIDE BY TWO AND RESTORE re_sign */
6644N/A+ z = z.divide_integer (2);
6644N/A+ return z.multiply_integer (re_sign);
6644N/A }
6644N/A
6644N/A /* Sets z = cosh x */
6644N/A public Number cosh ()
6644N/A {
6644N/A- var tmp = MPFloat.init2 ((Precision) precision);
6644N/A- tmp.cosh (re_num, Round.NEAREST);
6644N/A- var z = new Number.mpfloat (tmp);
6644N/A- return z;
6644N/A+ /* cosh (0) = 1 */
6644N/A+ if (is_zero ())
6644N/A+ return new Number.integer (1);
6644N/A+
6644N/A+ /* cosh (x) = (e^x + e^-x) / 2 */
6644N/A+ var t = abs ();
6644N/A+ t = t.epowy ();
6644N/A+ var z = t.reciprocal ();
6644N/A+ z = t.add (z);
6644N/A+ return z.divide_integer (2);
6644N/A }
6644N/A
6644N/A /* Sets z = tanh x */
6644N/A public Number tanh ()
6644N/A {
6644N/A- var tmp = MPFloat.init2 ((Precision) precision);
6644N/A- tmp.tanh (re_num, Round.NEAREST);
6644N/A- var z = new Number.mpfloat (tmp);
6644N/A+ /* tanh (0) = 0 */
6644N/A+ if (is_zero ())
6644N/A+ return new Number.integer (0);
6644N/A+
6644N/A+ var t = abs ();
6644N/A+
6644N/A+ /* SEE IF ABS (X) SO LARGE THAT RESULT IS +-1 */
6644N/A+ var r__1 = (float) T * 0.5 * Math.log ((float) BASE);
6644N/A+ var z = new Number.double (r__1);
6644N/A+ if (t.compare (z) > 0)
6644N/A+ return new Number.integer (re_sign);
6644N/A+
6644N/A+ /* If |x| >= 1/2 use ?, otherwise use ? to avoid cancellation */
6644N/A+ /* |tanh (x)| = (e^|2x| - 1) / (e^|2x| + 1) */
6644N/A+ t = t.multiply_integer (2);
6644N/A+ if (t.re_exponent > 0)
6644N/A+ {
6644N/A+ t = t.epowy ();
6644N/A+ z = t.add (new Number.integer (-1));
6644N/A+ t = t.add (new Number.integer (1));
6644N/A+ z = z.divide (t);
6644N/A+ }
6644N/A+ else
6644N/A+ {
6644N/A+ t = t.epowy ();
6644N/A+ z = t.add (new Number.integer (1));
6644N/A+ t = t.add (new Number.integer (-1));
6644N/A+ z = t.divide (z);
6644N/A+ }
6644N/A+
6644N/A+ /* Restore re_sign */
6644N/A+ z.re_sign = re_sign * z.re_sign;
6644N/A return z;
6644N/A }
6644N/A
6644N/A@@ -1014,10 +1496,12 @@
6644N/A /* Sets z = sinh⁻¹ x */
6644N/A public Number asinh ()
6644N/A {
6644N/A- var tmp = MPFloat.init2 ((Precision) precision);
6644N/A- tmp.asinh (re_num, Round.NEAREST);
6644N/A- var z = new Number.mpfloat (tmp);
6644N/A- return z;
6644N/A+ /* sinh⁻¹(x) = ln (x + √(x² + 1)) */
6644N/A+ var t = multiply (this);
6644N/A+ t = t.add (new Number.integer (1));
6644N/A+ t = t.sqrt ();
6644N/A+ t = add (t);
6644N/A+ return t.ln ();
6644N/A }
6644N/A
6644N/A /* Sets z = cosh⁻¹ x */
6644N/A@@ -1028,14 +1512,16 @@
6644N/A if (compare (t) < 0)
6644N/A {
6644N/A /* Translators: Error displayed when inverse hyperbolic cosine value is undefined */
6644N/A- error = _("Inverse hyperbolic cosine is undefined for values less than one");
6644N/A+ mperr (_("Inverse hyperbolic cosine is undefined for values less than one"));
6644N/A return new Number.integer (0);
6644N/A }
6644N/A
6644N/A- var tmp = MPFloat.init2 ((Precision) precision);
6644N/A- tmp.acosh (re_num, Round.NEAREST);
6644N/A- var z = new Number.mpfloat (tmp);
6644N/A- return z;
6644N/A+ /* cosh⁻¹(x) = ln (x + √(x² - 1)) */
6644N/A+ t = multiply (this);
6644N/A+ t = t.add (new Number.integer (-1));
6644N/A+ t = t.sqrt ();
6644N/A+ t = add (t);
6644N/A+ return t.ln ();
6644N/A }
6644N/A
6644N/A /* Sets z = tanh⁻¹ x */
6644N/A@@ -1045,14 +1531,17 @@
6644N/A if (compare (new Number.integer (1)) >= 0 || compare (new Number.integer (-1)) <= 0)
6644N/A {
6644N/A /* Translators: Error displayed when inverse hyperbolic tangent value is undefined */
6644N/A- error = _("Inverse hyperbolic tangent is undefined for values outside [-1, 1]");
6644N/A+ mperr (_("Inverse hyperbolic tangent is undefined for values outside [-1, 1]"));
6644N/A return new Number.integer (0);
6644N/A }
6644N/A
6644N/A- var tmp = MPFloat.init2 ((Precision) precision);
6644N/A- tmp.atanh (re_num, Round.NEAREST);
6644N/A- var z = new Number.mpfloat (tmp);
6644N/A- return z;
6644N/A+ /* atanh (x) = 0.5 * ln ((1 + x) / (1 - x)) */
6644N/A+ var n = add (new Number.integer (1));
6644N/A+ var d = invert_sign ();
6644N/A+ d = d.add (new Number.integer (1));
6644N/A+ var z = n.divide (d);
6644N/A+ z = z.ln ();
6644N/A+ return z.divide_integer (2);
6644N/A }
6644N/A
6644N/A /* Sets z = boolean AND for each bit in x and z */
6644N/A@@ -1062,7 +1551,7 @@
6644N/A is_positive_integer () || !y.is_positive_integer ())
6644N/A {
6644N/A /* Translators: Error displayed when boolean AND attempted on non-integer values */
6644N/A- error = _("Boolean AND is only defined for positive integers");
6644N/A+ mperr (_("Boolean AND is only defined for positive integers"));
6644N/A }
6644N/A
6644N/A return bitwise (y, (v1, v2) => { return v1 & v2; }, 0);
6644N/A@@ -1074,7 +1563,7 @@
6644N/A if (!is_positive_integer () || !y.is_positive_integer ())
6644N/A {
6644N/A /* Translators: Error displayed when boolean OR attempted on non-integer values */
6644N/A- error = _("Boolean OR is only defined for positive integers");
6644N/A+ mperr (_("Boolean OR is only defined for positive integers"));
6644N/A }
6644N/A
6644N/A return bitwise (y, (v1, v2) => { return v1 | v2; }, 0);
6644N/A@@ -1086,7 +1575,7 @@
6644N/A if (!is_positive_integer () || !y.is_positive_integer ())
6644N/A {
6644N/A /* Translators: Error displayed when boolean XOR attempted on non-integer values */
6644N/A- error = _("Boolean XOR is only defined for positive integers");
6644N/A+ mperr (_("Boolean XOR is only defined for positive integers"));
6644N/A }
6644N/A
6644N/A return bitwise (y, (v1, v2) => { return v1 ^ v2; }, 0);
6644N/A@@ -1098,7 +1587,7 @@
6644N/A if (!is_positive_integer ())
6644N/A {
6644N/A /* Translators: Error displayed when boolean XOR attempted on non-integer values */
6644N/A- error = _("Boolean NOT is only defined for positive integers");
6644N/A+ mperr (_("Boolean NOT is only defined for positive integers"));
6644N/A }
6644N/A
6644N/A return bitwise (new Number.integer (0), (v1, v2) => { return v1 ^ 0xF; }, wordlen);
6644N/A@@ -1121,7 +1610,7 @@
6644N/A if (!is_integer ())
6644N/A {
6644N/A /* Translators: Error displayed when bit shift attempted on non-integer values */
6644N/A- error = _("Shift is only possible on integer values");
6644N/A+ mperr (_("Shift is only possible on integer values"));
6644N/A return new Number.integer (0);
6644N/A }
6644N/A
6644N/A@@ -1253,60 +1742,1056 @@
6644N/A private Number copy ()
6644N/A {
6644N/A var z = new Number ();
6644N/A- var tmp = MPFloat.init2 ((Precision) precision);
6644N/A- var tmp2 = MPFloat.init2 ((Precision) precision);
6644N/A- tmp.@set(re_num, Round.NEAREST);
6644N/A- tmp2.@set(im_num, Round.NEAREST);
6644N/A- z.re_num = tmp;
6644N/A- z.im_num = tmp2;
6644N/A+ z.re_sign = re_sign;
6644N/A+ z.im_sign = im_sign;
6644N/A+ z.re_exponent = re_exponent;
6644N/A+ z.im_exponent = im_exponent;
6644N/A+ for (var i = 0; i < re_fraction.length; i++)
6644N/A+ {
6644N/A+ z.re_fraction[i] = re_fraction[i];
6644N/A+ z.im_fraction[i] = im_fraction[i];
6644N/A+ }
6644N/A+
6644N/A return z;
6644N/A }
6644N/A
6644N/A+ private Number add_with_sign (int y_sign, Number y)
6644N/A+ {
6644N/A+ if (is_complex () || y.is_complex ())
6644N/A+ {
6644N/A+ Number real_z, im_z;
6644N/A+
6644N/A+ var real_x = real_component ();
6644N/A+ var im_x = imaginary_component ();
6644N/A+ var real_y = y.real_component ();
6644N/A+ var im_y = y.imaginary_component ();
6644N/A+
6644N/A+ real_z = real_x.add_real (y_sign * y.re_sign, real_y);
6644N/A+ im_z = im_x.add_real (y_sign * y.im_sign, im_y);
6644N/A+
6644N/A+ return new Number.complex (real_z, im_z);
6644N/A+ }
6644N/A+ else
6644N/A+ return add_real (y_sign * y.re_sign, y);
6644N/A+ }
6644N/A+
6644N/A+
6644N/A private Number epowy_real ()
6644N/A {
6644N/A- var z = copy ();
6644N/A- var tmp = z.re_num;
6644N/A- tmp.exp (re_num, Round.NEAREST);
6644N/A- z.re_num = tmp;
6644N/A+ /* e^0 = 1 */
6644N/A+ if (is_zero ())
6644N/A+ return new Number.integer (1);
6644N/A+
6644N/A+ /* If |x| < 1 use exp */
6644N/A+ if (re_exponent <= 0)
6644N/A+ return exp ();
6644N/A+
6644N/A+ /* NOW SAFE TO CONVERT X TO REAL */
6644N/A+ var rx = to_double ();
6644N/A+
6644N/A+ /* SAVE re_sign AND WORK WITH ABS (X) */
6644N/A+ var xs = re_sign;
6644N/A+ var t2 = abs ();
6644N/A+
6644N/A+ /* GET fractional AND INTEGER PARTS OF ABS (X) */
6644N/A+ var ix = t2.to_integer ();
6644N/A+ t2 = t2.fractional_component ();
6644N/A+
6644N/A+ /* ATTACH re_sign TO fractional PART AND COMPUTE EXP OF IT */
6644N/A+ t2.re_sign *= xs;
6644N/A+ var z = t2.exp ();
6644N/A+
6644N/A+ /* COMPUTE E-2 OR 1/E USING TWO EXTRA DIGITS IN CASE ABS (X) LARGE
6644N/A+ * (BUT ONLY ONE EXTRA DIGIT IF T < 4)
6644N/A+ */
6644N/A+ var tss = 0;
6644N/A+ if (T < 4)
6644N/A+ tss = T + 1;
6644N/A+ else
6644N/A+ tss = T + 2;
6644N/A+
6644N/A+ /* LOOP FOR E COMPUTATION. DECREASE T IF POSSIBLE. */
6644N/A+ /* Computing MIN */
6644N/A+ var t1 = new Number.integer (xs);
6644N/A+
6644N/A+ t2.re_sign = 0;
6644N/A+ for (var i = 2 ; ; i++)
6644N/A+ {
6644N/A+ if (int.min (tss, tss + 2 + t1.re_exponent) <= 2)
6644N/A+ break;
6644N/A+
6644N/A+ t1 = t1.divide_integer (i * xs);
6644N/A+ t2 = t2.add (t1);
6644N/A+ if (t1.is_zero ())
6644N/A+ break;
6644N/A+ }
6644N/A+
6644N/A+ /* RAISE E OR 1/E TO POWER IX */
6644N/A+ if (xs > 0)
6644N/A+ t2 = t2.add (new Number.integer (2));
6644N/A+ t2 = t2.xpowy_integer (ix);
6644N/A+
6644N/A+ /* MULTIPLY EXPS OF INTEGER AND fractional PARTS */
6644N/A+ z = z.multiply (t2);
6644N/A+
6644N/A+ /* CHECK THAT RELATIVE ERROR LESS THAN 0.01 UNLESS ABS (X) LARGE
6644N/A+ * (WHEN EXP MIGHT OVERFLOW OR UNDERFLOW)
6644N/A+ */
6644N/A+ if (Math.fabs (rx) > 10.0f)
6644N/A+ return z;
6644N/A+ var rz = z.to_double ();
6644N/A+ var r__1 = rz - Math.exp (rx);
6644N/A+ if (Math.fabs (r__1) < rz * 0.01f)
6644N/A+ return z;
6644N/A+
6644N/A+ /* THE FOLLOWING MESSAGE MAY INDICATE THAT
6644N/A+ * B**(T-1) IS TOO SMALL, OR THAT M IS TOO SMALL SO THE
6644N/A+ * RESULT UNDERFLOWED.
6644N/A+ */
6644N/A+ mperr ("*** ERROR OCCURRED IN EPOWY, RESULT INCORRECT ***");
6644N/A return z;
6644N/A+
6644N/A }
6644N/A
6644N/A+ /* Return e^x for |x| < 1 USING AN o (SQRt (T).m (T)) ALGORITHM
6644N/A+ * DESCRIBED IN - R. P. BRENT, THE COMPLEXITY OF MULTIPLE-
6644N/A+ * PRECISION ARITHMETIC (IN COMPLEXITY OF COMPUTATIONAL PROBLEM
6644N/A+ * SOLVING, UNIV. OF QUEENSLAND PRESS, BRISBANE, 1976, 126-165).
6644N/A+ * ASYMPTOTICALLY FASTER METHODS EXIST, BUT ARE NOT USEFUL
6644N/A+ * UNLESS T IS VERY LARGE. SEE COMMENTS TO MP_ATAN AND MPPIGL.
6644N/A+ */
6644N/A+ private Number exp ()
6644N/A+ {
6644N/A+ /* e^0 = 1 */
6644N/A+ if (is_zero ())
6644N/A+ return new Number.integer (1);
6644N/A+
6644N/A+ /* Only defined for |x| < 1 */
6644N/A+ if (re_exponent > 0)
6644N/A+ {
6644N/A+ mperr ("*** ABS (X) NOT LESS THAN 1 IN CALL TO MP_EXP ***");
6644N/A+ return new Number.integer (0);
6644N/A+ }
6644N/A+
6644N/A+ var t1 = this;
6644N/A+ var rlb = Math.log (BASE);
6644N/A+
6644N/A+ /* Compute approximately optimal q (and divide x by 2^q) */
6644N/A+ var q = (int) (Math.sqrt (T * 0.48 * rlb) + re_exponent * 1.44 * rlb);
6644N/A+
6644N/A+ /* HALVE Q TIMES */
6644N/A+ if (q > 0)
6644N/A+ {
6644N/A+ var ib = BASE << 2;
6644N/A+ var ic = 1;
6644N/A+ for (var i = 1; i <= q; i++)
6644N/A+ {
6644N/A+ ic *= 2;
6644N/A+ if (ic < ib && ic != BASE && i < q)
6644N/A+ continue;
6644N/A+ t1 = t1.divide_integer (ic);
6644N/A+ ic = 1;
6644N/A+ }
6644N/A+ }
6644N/A+
6644N/A+ if (t1.is_zero ())
6644N/A+ return new Number.integer (0);
6644N/A+
6644N/A+ /* Sum series, reducing t where possible */
6644N/A+ var z = t1.copy ();
6644N/A+ var t2 = t1;
6644N/A+ for (var i = 2; T + t2.re_exponent - z.re_exponent > 0; i++)
6644N/A+ {
6644N/A+ t2 = t1.multiply (t2);
6644N/A+ t2 = t2.divide_integer (i);
6644N/A+ z = t2.add (z);
6644N/A+ if (t2.is_zero ())
6644N/A+ break;
6644N/A+ }
6644N/A+
6644N/A+ /* Apply (x+1)^2 - 1 = x (2 + x) for q iterations */
6644N/A+ for (var i = 1; i <= q; i++)
6644N/A+ {
6644N/A+ t1 = z.add (new Number.integer (2));
6644N/A+ z = t1.multiply (z);
6644N/A+ }
6644N/A+
6644N/A+ return z.add (new Number.integer (1));
6644N/A+ }
6644N/A+
6644N/A+ private Number pwr (Number y)
6644N/A+ {
6644N/A+ /* (-x)^y imaginary */
6644N/A+ /* FIXME: Make complex numbers optional */
6644N/A+ /*if (re_sign < 0)
6644N/A+ {
6644N/A+ mperr (_("The power of negative numbers is only defined for integer exponents"));
6644N/A+ return new Number.integer (0);
6644N/A+ }*/
6644N/A+
6644N/A+ /* 0^y = 0, 0^-y undefined */
6644N/A+ if (is_zero ())
6644N/A+ {
6644N/A+ if (y.re_sign < 0)
6644N/A+ mperr (_("The power of zero is undefined for a negative exponent"));
6644N/A+ return new Number.integer (0);
6644N/A+ }
6644N/A+
6644N/A+ /* x^0 = 1 */
6644N/A+ if (y.is_zero ())
6644N/A+ return new Number.integer (1);
6644N/A+
6644N/A+ return y.multiply (ln ()).epowy ();
6644N/A+ }
6644N/A+
6644N/A private Number root_real (int64 n)
6644N/A {
6644N/A+ /* x^(1/1) = x */
6644N/A+ if (n == 1)
6644N/A+ return this;
6644N/A+
6644N/A+ /* x^(1/0) invalid */
6644N/A+ if (n == 0)
6644N/A+ {
6644N/A+ mperr (_("Root must be non-zero"));
6644N/A+ return new Number.integer (0);
6644N/A+ }
6644N/A+
6644N/A+ var np = n.abs ();
6644N/A+
6644N/A+ /* LOSS OF ACCURACY IF NP LARGE, SO ONLY ALLOW NP <= MAX (B, 64) */
6644N/A+ if (np > int.max (BASE, 64))
6644N/A+ {
6644N/A+ mperr ("*** ABS (N) TOO LARGE IN CALL TO ROOT ***");
6644N/A+ return new Number.integer (0);
6644N/A+ }
6644N/A+
6644N/A+ /* 0^(1/n) = 0 for positive n */
6644N/A+ if (is_zero ())
6644N/A+ {
6644N/A+ if (n <= 0)
6644N/A+ mperr (_("Negative root of zero is undefined"));
6644N/A+ return new Number.integer (0);
6644N/A+ }
6644N/A+
6644N/A+ // FIXME: Imaginary root
6644N/A+ if (re_sign < 0 && np % 2 == 0)
6644N/A+ {
6644N/A+ mperr (_("nth root of negative number is undefined for even n"));
6644N/A+ return new Number.integer (0);
6644N/A+ }
6644N/A+
6644N/A+ /* DIVIDE re_exponent BY NP */
6644N/A+ var ex = re_exponent / np;
6644N/A+
6644N/A+ /* Get initial approximation */
6644N/A+ var t1 = copy ();
6644N/A+ t1.re_exponent = 0;
6644N/A+ var approximation = Math.exp (((np * ex - re_exponent) * Math.log (BASE) - Math.log (Math.fabs (t1.to_float ()))) / np);
6644N/A+ t1 = new Number.double (approximation);
6644N/A+ t1.re_sign = re_sign;
6644N/A+ t1.re_exponent -= (int) ex;
6644N/A+
6644N/A+ /* MAIN ITERATION LOOP */
6644N/A+ var it0 = 3;
6644N/A+ var t = it0;
6644N/A+ Number t2;
6644N/A+ while (true)
6644N/A+ {
6644N/A+ /* t1 = t1 - ((t1 * ((x * t1^np) - 1)) / np) */
6644N/A+ t2 = t1.xpowy_integer (np);
6644N/A+ t2 = multiply (t2);
6644N/A+ t2 = t2.add (new Number.integer (-1));
6644N/A+ t2 = t1.multiply (t2);
6644N/A+ t2 = t2.divide_integer (np);
6644N/A+ t1 = t1.subtract (t2);
6644N/A+
6644N/A+ /* FOLLOWING LOOP ALMOST DOUBLES T (POSSIBLE BECAUSE
6644N/A+ * NEWTONS METHOD HAS 2ND ORDER CONVERGENCE).
6644N/A+ */
6644N/A+ if (t >= T)
6644N/A+ break;
6644N/A+
6644N/A+ var ts3 = t;
6644N/A+ var ts2 = t;
6644N/A+ t = T;
6644N/A+ do
6644N/A+ {
6644N/A+ ts2 = t;
6644N/A+ t = (t + it0) / 2;
6644N/A+ } while (t > ts3);
6644N/A+ t = int.min (ts2, T);
6644N/A+ }
6644N/A+
6644N/A+ /* NOW r (I2) IS X**(-1/NP)
6644N/A+ * CHECK THAT NEWTON ITERATION WAS CONVERGING
6644N/A+ */
6644N/A+ if (t2.re_sign != 0 && (t1.re_exponent - t2.re_exponent) << 1 < T - it0)
6644N/A+ {
6644N/A+ /* THE FOLLOWING MESSAGE MAY INDICATE THAT B**(T-1) IS TOO SMALL,
6644N/A+ * OR THAT THE INITIAL APPROXIMATION OBTAINED USING ALOG AND EXP
6644N/A+ * IS NOT ACCURATE ENOUGH.
6644N/A+ */
6644N/A+ mperr ("*** ERROR OCCURRED IN ROOT, NEWTON ITERATION NOT CONVERGING PROPERLY ***");
6644N/A+ }
6644N/A+
6644N/A+ if (n >= 0)
6644N/A+ {
6644N/A+ t1 = t1.xpowy_integer (n - 1);
6644N/A+ return multiply (t1);
6644N/A+ }
6644N/A+
6644N/A+ return t1;
6644N/A+
6644N/A+ }
6644N/A+
6644N/A+ /* ROUTINE CALLED BY DIVIDE AND SQRT TO ENSURE THAT
6644N/A+ * RESULTS ARE REPRESENTED EXACTLY IN T-2 DIGITS IF THEY
6644N/A+ * CAN BE. X IS AN MP NUMBER, I AND J ARE INTEGERS.
6644N/A+ */
6644N/A+ private Number ext (int i, int j)
6644N/A+ {
6644N/A+ if (is_zero () || T <= 2 || i == 0)
6644N/A+ return this;
6644N/A+
6644N/A+ /* COMPUTE MAXIMUM POSSIBLE ERROR IN THE LAST PLACE */
6644N/A+ var q = (j + 1) / i + 1;
6644N/A+ var s = BASE * re_fraction[T - 2] + re_fraction[T - 1];
6644N/A+
6644N/A+ /* SET LAST TWO DIGITS TO ZERO */
6644N/A var z = copy ();
6644N/A- var tmp = z.re_num;
6644N/A- tmp.root (re_num, (ulong) n, Round.NEAREST);
6644N/A- z.re_num = tmp;
6644N/A- return z;
6644N/A+ if (s <= q)
6644N/A+ {
6644N/A+ z.re_fraction[T - 2] = 0;
6644N/A+ z.re_fraction[T - 1] = 0;
6644N/A+ return z;
6644N/A+ }
6644N/A+
6644N/A+ if (s + q < BASE * BASE)
6644N/A+ return z;
6644N/A+
6644N/A+ /* ROUND UP HERE */
6644N/A+ z.re_fraction[T - 2] = BASE - 1;
6644N/A+ z.re_fraction[T - 1] = BASE;
6644N/A+
6644N/A+ /* NORMALIZE X (LAST DIGIT B IS OK IN MULTIPLY_INTEGER) */
6644N/A+ return z.multiply_integer (1);
6644N/A }
6644N/A
6644N/A+
6644N/A private Number ln_real ()
6644N/A {
6644N/A- var tmp = MPFloat.init2 ((Precision) precision);
6644N/A- tmp.log (re_num, Round.NEAREST);
6644N/A- var z = new Number.mpfloat (tmp);
6644N/A+ // ln(e^1) = 1, fixes precision loss
6644N/A+ if (equals (new Number.eulers ()))
6644N/A+ return new Number.integer (1);
6644N/A+
6644N/A+ /* LOOP TO GET APPROXIMATE Ln (X) USING SINGLE-PRECISION */
6644N/A+ var t1 = copy ();
6644N/A+ var z = new Number.integer (0);
6644N/A+ for (var k = 0; k < 10; k++)
6644N/A+ {
6644N/A+ /* COMPUTE FINAL CORRECTION ACCURATELY USING LNS */
6644N/A+ var t2 = t1.add (new Number.integer (-1));
6644N/A+ if (t2.is_zero () || t2.re_exponent + 1 <= 0)
6644N/A+ return z.add (t2.lns ());
6644N/A+
6644N/A+ /* REMOVE EXPONENT TO AVOID FLOATING-POINT OVERFLOW */
6644N/A+ var e = t1.re_exponent;
6644N/A+ t1.re_exponent = 0;
6644N/A+ var rx = t1.to_double ();
6644N/A+ t1.re_exponent = e;
6644N/A+ var rlx = Math.log (rx) + e * Math.log (BASE);
6644N/A+ t2 = new Number.double (-rlx);
6644N/A+
6644N/A+ /* UPDATE Z AND COMPUTE ACCURATE EXP OF APPROXIMATE LOG */
6644N/A+ z = z.subtract (t2);
6644N/A+ t2 = t2.epowy ();
6644N/A+
6644N/A+ /* COMPUTE RESIDUAL WHOSE LOG IS STILL TO BE FOUND */
6644N/A+ t1 = t1.multiply (t2);
6644N/A+ }
6644N/A+
6644N/A+ mperr ("*** ERROR IN LN, ITERATION NOT CONVERGING ***");
6644N/A return z;
6644N/A+
6644N/A }
6644N/A
6644N/A+ /* RETURNS MP Y = Ln (1+X) IF X IS AN MP NUMBER SATISFYING THE
6644N/A+ * CONDITION ABS (X) < 1/B, ERROR OTHERWISE.
6644N/A+ * USES NEWTONS METHOD TO SOLVE THE EQUATION
6644N/A+ * EXP1(-Y) = X, THEN REVERSES re_sign OF Y.
6644N/A+ */
6644N/A+ private Number lns ()
6644N/A+ {
6644N/A+ /* ln (1+0) = 0 */
6644N/A+ if (is_zero ())
6644N/A+ return this;
6644N/A+
6644N/A+ /* Get starting approximation -ln (1+x) ~= -x + x^2/2 - x^3/3 + x^4/4 */
6644N/A+ var t2 = copy ();
6644N/A+ var t1 = divide_integer (4);
6644N/A+ t1 = t1.add (new Number.fraction (-1, 3));
6644N/A+ t1 = multiply (t1);
6644N/A+ t1 = t1.add (new Number.fraction (1, 2));
6644N/A+ t1 = multiply (t1);
6644N/A+ t1 = t1.add (new Number.integer (-1));
6644N/A+ var z = multiply (t1);
6644N/A+
6644N/A+ /* Solve using Newtons method */
6644N/A+ var it0 = 5;
6644N/A+ var t = it0;
6644N/A+ Number t3;
6644N/A+ while (true)
6644N/A+ {
6644N/A+ /* t3 = (e^t3 - 1) */
6644N/A+ /* z = z - (t2 + t3 + (t2 * t3)) */
6644N/A+ t3 = z.epowy ();
6644N/A+ t3 = t3.add (new Number.integer (-1));
6644N/A+ t1 = t2.multiply (t3);
6644N/A+ t3 = t3.add (t1);
6644N/A+ t3 = t2.add (t3);
6644N/A+ z = z.subtract (t3);
6644N/A+ if (t >= T)
6644N/A+ break;
6644N/A+
6644N/A+ /* FOLLOWING LOOP COMPUTES NEXT VALUE OF T TO USE.
6644N/A+ * BECAUSE NEWTONS METHOD HAS 2ND ORDER CONVERGENCE,
6644N/A+ * WE CAN ALMOST DOUBLE T EACH TIME.
6644N/A+ */
6644N/A+ var ts3 = t;
6644N/A+ var ts2 = t;
6644N/A+ t = T;
6644N/A+ do
6644N/A+ {
6644N/A+ ts2 = t;
6644N/A+ t = (t + it0) / 2;
6644N/A+ } while (t > ts3);
6644N/A+ t = ts2;
6644N/A+ }
6644N/A+
6644N/A+ /* CHECK THAT NEWTON ITERATION WAS CONVERGING AS EXPECTED */
6644N/A+ if (t3.re_sign != 0 && t3.re_exponent << 1 > it0 - T)
6644N/A+ mperr ("*** ERROR OCCURRED IN LNS, NEWTON ITERATION NOT CONVERGING PROPERLY ***");
6644N/A+
6644N/A+ z.re_sign = -z.re_sign;
6644N/A+
6644N/A+ return z;
6644N/A+ }
6644N/A+
6644N/A+ private Number add_real (int y_sign, Number y)
6644N/A+ {
6644N/A+ var re_sign_prod = y_sign * re_sign;
6644N/A+
6644N/A+ /* 0 + y = y */
6644N/A+ if (is_zero ())
6644N/A+ {
6644N/A+ if (y_sign != y.re_sign)
6644N/A+ return y.invert_sign ();
6644N/A+ else
6644N/A+ return y;
6644N/A+ }
6644N/A+ /* x + 0 = x */
6644N/A+ else if (y.is_zero ())
6644N/A+ return this;
6644N/A+
6644N/A+ var exp_diff = re_exponent - y.re_exponent;
6644N/A+ var med = exp_diff.abs ();
6644N/A+ var x_largest = false;
6644N/A+ if (exp_diff < 0)
6644N/A+ x_largest = false;
6644N/A+ else if (exp_diff > 0)
6644N/A+ x_largest = true;
6644N/A+ else
6644N/A+ {
6644N/A+ /* EXPONENTS EQUAL SO COMPARE SIGNS, THEN FRACTIONS IF NEC. */
6644N/A+ if (re_sign_prod < 0)
6644N/A+ {
6644N/A+ /* signs are not equal. find out which mantissa is larger. */
6644N/A+ int j;
6644N/A+ for (j = 0; j < T; j++)
6644N/A+ {
6644N/A+ int i = re_fraction[j] - y.re_fraction[j];
6644N/A+ if (i != 0)
6644N/A+ {
6644N/A+ if (i < 0)
6644N/A+ x_largest = false;
6644N/A+ else if (i > 0)
6644N/A+ x_largest = true;
6644N/A+ break;
6644N/A+ }
6644N/A+ }
6644N/A+
6644N/A+ /* Both mantissas equal, so result is zero. */
6644N/A+ if (j >= T)
6644N/A+ return new Number.integer (0);
6644N/A+ }
6644N/A+ }
6644N/A+
6644N/A+ var z = new Number.integer (0);
6644N/A+ int[] big_fraction, small_fraction;
6644N/A+ if (x_largest)
6644N/A+ {
6644N/A+ z.re_sign = re_sign;
6644N/A+ z.re_exponent = re_exponent;
6644N/A+ big_fraction = re_fraction;
6644N/A+ small_fraction = y.re_fraction;
6644N/A+ }
6644N/A+ else
6644N/A+ {
6644N/A+ z.re_sign = y_sign;
6644N/A+ z.re_exponent = y.re_exponent;
6644N/A+ big_fraction = y.re_fraction;
6644N/A+ small_fraction = re_fraction;
6644N/A+ }
6644N/A+
6644N/A+ /* CLEAR GUARD DIGITS TO RIGHT OF X DIGITS */
6644N/A+ for (var i = 3; i >= med; i--)
6644N/A+ z.re_fraction[T + i] = 0;
6644N/A+
6644N/A+ if (re_sign_prod >= 0)
6644N/A+ {
6644N/A+ /* This is probably insufficient overflow detection, but it makes us not crash at least. */
6644N/A+ if (T + 3 < med)
6644N/A+ {
6644N/A+ mperr (_("Overflow: the result couldn't be calculated"));
6644N/A+ return new Number.integer (0);
6644N/A+ }
6644N/A+
6644N/A+ /* HERE DO ADDITION, re_exponent (Y) >= re_exponent (X) */
6644N/A+ var i = 0;
6644N/A+ for (i = T + 3; i >= T; i--)
6644N/A+ z.re_fraction[i] = small_fraction[i - med];
6644N/A+
6644N/A+ var c = 0;
6644N/A+ for (; i >= med; i--)
6644N/A+ {
6644N/A+ c = big_fraction[i] + small_fraction[i - med] + c;
6644N/A+
6644N/A+ if (c < BASE)
6644N/A+ {
6644N/A+ /* NO CARRY GENERATED HERE */
6644N/A+ z.re_fraction[i] = c;
6644N/A+ c = 0;
6644N/A+ }
6644N/A+ else
6644N/A+ {
6644N/A+ /* CARRY GENERATED HERE */
6644N/A+ z.re_fraction[i] = c - BASE;
6644N/A+ c = 1;
6644N/A+ }
6644N/A+ }
6644N/A+
6644N/A+ for (; i >= 0; i--)
6644N/A+ {
6644N/A+ c = big_fraction[i] + c;
6644N/A+ if (c < BASE)
6644N/A+ {
6644N/A+ z.re_fraction[i] = c;
6644N/A+ i--;
6644N/A+
6644N/A+ /* NO CARRY POSSIBLE HERE */
6644N/A+ for (; i >= 0; i--)
6644N/A+ z.re_fraction[i] = big_fraction[i];
6644N/A+
6644N/A+ c = 0;
6644N/A+ break;
6644N/A+ }
6644N/A+
6644N/A+ z.re_fraction[i] = 0;
6644N/A+ c = 1;
6644N/A+ }
6644N/A+
6644N/A+ /* MUST SHIFT RIGHT HERE AS CARRY OFF END */
6644N/A+ if (c != 0)
6644N/A+ {
6644N/A+ for (var j = T + 3; j > 0; j--)
6644N/A+ z.re_fraction[j] = z.re_fraction[j - 1];
6644N/A+ z.re_fraction[0] = 1;
6644N/A+ z.re_exponent++;
6644N/A+ }
6644N/A+ }
6644N/A+ else
6644N/A+ {
6644N/A+ var c = 0;
6644N/A+ var i = 0;
6644N/A+ for (i = T + med - 1; i >= T; i--)
6644N/A+ {
6644N/A+ /* HERE DO SUBTRACTION, ABS (Y) > ABS (X) */
6644N/A+ z.re_fraction[i] = c - small_fraction[i - med];
6644N/A+ c = 0;
6644N/A+
6644N/A+ /* BORROW GENERATED HERE */
6644N/A+ if (z.re_fraction[i] < 0)
6644N/A+ {
6644N/A+ c = -1;
6644N/A+ z.re_fraction[i] += BASE;
6644N/A+ }
6644N/A+ }
6644N/A+
6644N/A+ for (; i >= med; i--)
6644N/A+ {
6644N/A+ c = big_fraction[i] + c - small_fraction[i - med];
6644N/A+ if (c >= 0)
6644N/A+ {
6644N/A+ /* NO BORROW GENERATED HERE */
6644N/A+ z.re_fraction[i] = c;
6644N/A+ c = 0;
6644N/A+ }
6644N/A+ else
6644N/A+ {
6644N/A+ /* BORROW GENERATED HERE */
6644N/A+ z.re_fraction[i] = c + BASE;
6644N/A+ c = -1;
6644N/A+ }
6644N/A+ }
6644N/A+
6644N/A+ for (; i >= 0; i--)
6644N/A+ {
6644N/A+ c = big_fraction[i] + c;
6644N/A+
6644N/A+ if (c >= 0)
6644N/A+ {
6644N/A+ z.re_fraction[i] = c;
6644N/A+ i--;
6644N/A+
6644N/A+ /* NO CARRY POSSIBLE HERE */
6644N/A+ for (; i >= 0; i--)
6644N/A+ z.re_fraction[i] = big_fraction[i];
6644N/A+
6644N/A+ break;
6644N/A+ }
6644N/A+
6644N/A+ z.re_fraction[i] = c + BASE;
6644N/A+ c = -1;
6644N/A+ }
6644N/A+ }
6644N/A+
6644N/A+ mp_normalize (ref z);
6644N/A+
6644N/A+ return z;
6644N/A+ }
6644N/A+
6644N/A+ private Number multiply_real (Number y)
6644N/A+ {
6644N/A+ /* x*0 = 0*y = 0 */
6644N/A+ if (re_sign == 0 || y.re_sign == 0)
6644N/A+ return new Number.integer (0);
6644N/A+
6644N/A+ var z = new Number.integer (0);
6644N/A+ z.re_sign = re_sign * y.re_sign;
6644N/A+ z.re_exponent = re_exponent + y.re_exponent;
6644N/A+
6644N/A+ var r = new Number.integer (0);
6644N/A+
6644N/A+ /* PERFORM MULTIPLICATION */
6644N/A+ var c = 8;
6644N/A+ for (var i = 0; i < T; i++)
6644N/A+ {
6644N/A+ var xi = re_fraction[i];
6644N/A+
6644N/A+ /* FOR SPEED, PUT THE NUMBER WITH MANY ZEROS FIRST */
6644N/A+ if (xi == 0)
6644N/A+ continue;
6644N/A+
6644N/A+ /* Computing MIN */
6644N/A+ for (var j = 0; j < int.min (T, T + 3 - i); j++)
6644N/A+ r.re_fraction[i+j+1] += xi * y.re_fraction[j];
6644N/A+ c--;
6644N/A+ if (c > 0)
6644N/A+ continue;
6644N/A+
6644N/A+ /* CHECK FOR LEGAL BASE B DIGIT */
6644N/A+ if (xi < 0 || xi >= BASE)
6644N/A+ {
6644N/A+ mperr ("*** ILLEGAL BASE B DIGIT IN CALL TO MULTIPLY, POSSIBLE OVERWRITING PROBLEM ***");
6644N/A+ return new Number.integer (0);
6644N/A+ }
6644N/A+
6644N/A+ /* PROPAGATE CARRIES AT END AND EVERY EIGHTH TIME,
6644N/A+ * FASTER THAN DOING IT EVERY TIME.
6644N/A+ */
6644N/A+ for (var j = T + 3; j >= 0; j--)
6644N/A+ {
6644N/A+ int ri = r.re_fraction[j] + c;
6644N/A+ if (ri < 0)
6644N/A+ {
6644N/A+ mperr ("*** INTEGER OVERFLOW IN MULTIPLY, B TOO LARGE ***");
6644N/A+ return new Number.integer (0);
6644N/A+ }
6644N/A+ c = ri / BASE;
6644N/A+ r.re_fraction[j] = ri - BASE * c;
6644N/A+ }
6644N/A+ if (c != 0)
6644N/A+ {
6644N/A+ mperr ("*** ILLEGAL BASE B DIGIT IN CALL TO MULTIPLY, POSSIBLE OVERWRITING PROBLEM ***");
6644N/A+ return new Number.integer (0);
6644N/A+ }
6644N/A+ c = 8;
6644N/A+ }
6644N/A+
6644N/A+ if (c != 8)
6644N/A+ {
6644N/A+ c = 0;
6644N/A+ for (var i = T + 3; i >= 0; i--)
6644N/A+ {
6644N/A+ int ri = r.re_fraction[i] + c;
6644N/A+ if (ri < 0)
6644N/A+ {
6644N/A+ mperr ("*** INTEGER OVERFLOW IN MULTIPLY, B TOO LARGE ***");
6644N/A+ return new Number.integer (0);
6644N/A+ }
6644N/A+ c = ri / BASE;
6644N/A+ r.re_fraction[i] = ri - BASE * c;
6644N/A+ }
6644N/A+
6644N/A+ if (c != 0)
6644N/A+ {
6644N/A+ mperr ("*** ILLEGAL BASE B DIGIT IN CALL TO MULTIPLY, POSSIBLE OVERWRITING PROBLEM ***");
6644N/A+ return new Number.integer (0);
6644N/A+ }
6644N/A+ }
6644N/A+
6644N/A+ /* Clear complex part */
6644N/A+ z.im_sign = 0;
6644N/A+ z.im_exponent = 0;
6644N/A+ for (var i = 0; i < z.im_fraction.length; i++)
6644N/A+ z.im_fraction[i] = 0;
6644N/A+
6644N/A+ /* NORMALIZE AND ROUND RESULT */
6644N/A+ // FIXME: Use stack variable because of mp_normalize brokeness
6644N/A+ for (var i = 0; i < SIZE; i++)
6644N/A+ z.re_fraction[i] = r.re_fraction[i];
6644N/A+ mp_normalize (ref z);
6644N/A+
6644N/A+ return z;
6644N/A+ }
6644N/A+
6644N/A+ private Number multiply_integer_real (int64 y)
6644N/A+ {
6644N/A+ /* x*0 = 0*y = 0 */
6644N/A+ if (is_zero () || y == 0)
6644N/A+ return new Number.integer (0);
6644N/A+
6644N/A+ /* x*1 = x, x*-1 = -x */
6644N/A+ // FIXME: Why is this not working? mp_ext is using this function to do a normalization
6644N/A+ /*if (y == 1 || y == -1)
6644N/A+ {
6644N/A+ if (y < 0)
6644N/A+ z = invert_sign ();
6644N/A+ else
6644N/A+ z = this;
6644N/A+ return z;
6644N/A+ }*/
6644N/A+
6644N/A+ /* Copy x as z may also refer to x */
6644N/A+ var z = new Number.integer (0);
6644N/A+
6644N/A+ if (y < 0)
6644N/A+ {
6644N/A+ y = -y;
6644N/A+ z.re_sign = -re_sign;
6644N/A+ }
6644N/A+ else
6644N/A+ z.re_sign = re_sign;
6644N/A+ z.re_exponent = re_exponent + 4;
6644N/A+
6644N/A+ /* FORM PRODUCT IN ACCUMULATOR */
6644N/A+ int64 c = 0;
6644N/A+
6644N/A+ /* IF y*B NOT REPRESENTABLE AS AN INTEGER WE HAVE TO SIMULATE
6644N/A+ * DOUBLE-PRECISION MULTIPLICATION.
6644N/A+ */
6644N/A+
6644N/A+ /* Computing MAX */
6644N/A+ if (y >= int.max (BASE << 3, 32767 / BASE))
6644N/A+ {
6644N/A+ /* HERE J IS TOO LARGE FOR SINGLE-PRECISION MULTIPLICATION */
6644N/A+ var j1 = y / BASE;
6644N/A+ var j2 = y - j1 * BASE;
6644N/A+
6644N/A+ /* FORM PRODUCT */
6644N/A+ for (var i = T + 3; i >= 0; i--)
6644N/A+ {
6644N/A+ var c1 = c / BASE;
6644N/A+ var c2 = c - BASE * c1;
6644N/A+ var ix = 0;
6644N/A+ if (i > 3)
6644N/A+ ix = re_fraction[i - 4];
6644N/A+
6644N/A+ var t = j2 * ix + c2;
6644N/A+ var is = t / BASE;
6644N/A+ c = j1 * ix + c1 + is;
6644N/A+ z.re_fraction[i] = (int) (t - BASE * is);
6644N/A+ }
6644N/A+ }
6644N/A+ else
6644N/A+ {
6644N/A+ int64 ri = 0;
6644N/A+ for (var i = T + 3; i >= 4; i--)
6644N/A+ {
6644N/A+ ri = y * re_fraction[i - 4] + c;
6644N/A+ c = ri / BASE;
6644N/A+ z.re_fraction[i] = (int) (ri - BASE * c);
6644N/A+ }
6644N/A+
6644N/A+ /* CHECK FOR INTEGER OVERFLOW */
6644N/A+ if (ri < 0)
6644N/A+ {
6644N/A+ mperr ("*** INTEGER OVERFLOW IN multiply_integer, B TOO LARGE ***");
6644N/A+ return new Number.integer (0);
6644N/A+ }
6644N/A+
6644N/A+ /* HAVE TO TREAT FIRST FOUR WORDS OF R SEPARATELY */
6644N/A+ for (var i = 3; i >= 0; i--)
6644N/A+ {
6644N/A+ var t = c;
6644N/A+ c = t / BASE;
6644N/A+ z.re_fraction[i] = (int) (t - BASE * c);
6644N/A+ }
6644N/A+ }
6644N/A+
6644N/A+ /* HAVE TO SHIFT RIGHT HERE AS CARRY OFF END */
6644N/A+ while (c != 0)
6644N/A+ {
6644N/A+ for (var i = T + 3; i >= 1; i--)
6644N/A+ z.re_fraction[i] = z.re_fraction[i - 1];
6644N/A+ var t = c;
6644N/A+ c = t / BASE;
6644N/A+ z.re_fraction[0] = (int) (t - BASE * c);
6644N/A+ z.re_exponent++;
6644N/A+ }
6644N/A+
6644N/A+ if (c < 0)
6644N/A+ {
6644N/A+ mperr ("*** INTEGER OVERFLOW IN multiply_integer, B TOO LARGE ***");
6644N/A+ return new Number.integer (0);
6644N/A+ }
6644N/A+
6644N/A+ z.im_sign = 0;
6644N/A+ z.im_exponent = 0;
6644N/A+ for (var i = 0; i < z.im_fraction.length; i++)
6644N/A+ z.im_fraction[i] = 0;
6644N/A+ mp_normalize (ref z);
6644N/A+ return z;
6644N/A+ }
6644N/A+
6644N/A private Number reciprocal_real ()
6644N/A {
6644N/A /* 1/0 invalid */
6644N/A if (is_zero ())
6644N/A {
6644N/A- error = _("Reciprocal of zero is undefined");
6644N/A+ mperr (_("Reciprocal of zero is undefined"));
6644N/A return new Number.integer (0);
6644N/A }
6644N/A
6644N/A- var tmp = MPFloat.init2 ((Precision) precision);
6644N/A- tmp.set_unsigned_integer (1, Round.NEAREST);
6644N/A- tmp.divide (tmp, re_num, Round.NEAREST);
6644N/A- var z = copy ();
6644N/A- z.re_num.clear ();
6644N/A- z.re_num = tmp;
6644N/A- return z;
6644N/A+ /* Start by approximating value using floating point */
6644N/A+ var t1 = copy ();
6644N/A+ t1.re_exponent = 0;
6644N/A+ t1 = new Number.double (1.0 / t1.to_double ());
6644N/A+ t1.re_exponent -= re_exponent;
6644N/A+
6644N/A+ var t = 3;
6644N/A+ var it0 = t;
6644N/A+ Number t2;
6644N/A+ while (true)
6644N/A+ {
6644N/A+ /* t1 = t1 - (t1 * ((x * t1) - 1)) (2*t1 - t1^2*x) */
6644N/A+ t2 = multiply (t1);
6644N/A+ t2 = t2.add (new Number.integer (-1));
6644N/A+ t2 = t1.multiply (t2);
6644N/A+ t1 = t1.subtract (t2);
6644N/A+ if (t >= T)
6644N/A+ break;
6644N/A+
6644N/A+ /* FOLLOWING LOOP ALMOST DOUBLES T (POSSIBLE
6644N/A+ * BECAUSE NEWTONS METHOD HAS 2ND ORDER CONVERGENCE).
6644N/A+ */
6644N/A+ var ts3 = t;
6644N/A+ var ts2 = 0;
6644N/A+ t = T;
6644N/A+ do
6644N/A+ {
6644N/A+ ts2 = t;
6644N/A+ t = (t + it0) / 2;
6644N/A+ } while (t > ts3);
6644N/A+ t = int.min (ts2, T);
6644N/A+ }
6644N/A+
6644N/A+ /* RETURN IF NEWTON ITERATION WAS CONVERGING */
6644N/A+ if (t2.re_sign != 0 && (t1.re_exponent - t2.re_exponent) << 1 < T - it0)
6644N/A+ {
6644N/A+ /* THE FOLLOWING MESSAGE MAY INDICATE THAT B**(T-1) IS TOO SMALL,
6644N/A+ * OR THAT THE STARTING APPROXIMATION IS NOT ACCURATE ENOUGH.
6644N/A+ */
6644N/A+ mperr ("*** ERROR OCCURRED IN RECIPROCAL, NEWTON ITERATION NOT CONVERGING PROPERLY ***");
6644N/A+ }
6644N/A+
6644N/A+ return t1;
6644N/A }
6644N/A
6644N/A+ private Number divide_integer_real (int64 y)
6644N/A+ {
6644N/A+ /* x/0 */
6644N/A+ if (y == 0)
6644N/A+ {
6644N/A+ /* Translators: Error displayed attempted to divide by zero */
6644N/A+ mperr (_("Division by zero is undefined"));
6644N/A+ return new Number.integer (0);
6644N/A+ }
6644N/A
6644N/A+ /* 0/y = 0 */
6644N/A+ if (is_zero ())
6644N/A+ return new Number.integer (0);
6644N/A+
6644N/A+ /* Division by -1 or 1 just changes re_sign */
6644N/A+ if (y == 1 || y == -1)
6644N/A+ {
6644N/A+ if (y < 0)
6644N/A+ return invert_sign ();
6644N/A+ else
6644N/A+ return this;
6644N/A+ }
6644N/A+
6644N/A+ var z = new Number.integer (0);
6644N/A+ if (y < 0)
6644N/A+ {
6644N/A+ y = -y;
6644N/A+ z.re_sign = -re_sign;
6644N/A+ }
6644N/A+ else
6644N/A+ z.re_sign = re_sign;
6644N/A+ z.re_exponent = re_exponent;
6644N/A+
6644N/A+ int64 c = 0;
6644N/A+ int64 i = 0;
6644N/A+
6644N/A+ /* IF y*B NOT REPRESENTABLE AS AN INTEGER HAVE TO SIMULATE
6644N/A+ * LONG DIVISION. ASSUME AT LEAST 16-BIT WORD.
6644N/A+ */
6644N/A+
6644N/A+ /* Computing MAX */
6644N/A+ var b2 = int.max (BASE << 3, 32767 / BASE);
6644N/A+ if (y < b2)
6644N/A+ {
6644N/A+ /* LOOK FOR FIRST NONZERO DIGIT IN QUOTIENT */
6644N/A+ int64 r1 = 0;
6644N/A+ do
6644N/A+ {
6644N/A+ c = BASE * c;
6644N/A+ if (i < T)
6644N/A+ c += re_fraction[i];
6644N/A+ i++;
6644N/A+ r1 = c / y;
6644N/A+ if (r1 < 0)
6644N/A+ {
6644N/A+ mperr ("*** INTEGER OVERFLOW IN DIVIDE_INTEGER, B TOO LARGE ***");
6644N/A+ return new Number.integer (0);
6644N/A+ }
6644N/A+ } while (r1 == 0);
6644N/A+
6644N/A+ /* ADJUST re_exponent AND GET T+4 DIGITS IN QUOTIENT */
6644N/A+ z.re_exponent += (int) (1 - i);
6644N/A+ z.re_fraction[0] = (int) r1;
6644N/A+ c = BASE * (c - y * r1);
6644N/A+ int64 kh = 1;
6644N/A+ if (i < T)
6644N/A+ {
6644N/A+ kh = T + 1 - i;
6644N/A+ for (var k = 1; k < kh; k++)
6644N/A+ {
6644N/A+ c += re_fraction[i];
6644N/A+ z.re_fraction[k] = (int) (c / y);
6644N/A+ c = BASE * (c - y * z.re_fraction[k]);
6644N/A+ i++;
6644N/A+ }
6644N/A+ if (c < 0)
6644N/A+ {
6644N/A+ mperr ("*** INTEGER OVERFLOW IN DIVIDE_INTEGER, B TOO LARGE ***");
6644N/A+ return new Number.integer (0);
6644N/A+ }
6644N/A+ }
6644N/A+
6644N/A+ for (var k = kh; k < T + 4; k++)
6644N/A+ {
6644N/A+ z.re_fraction[k] = (int) (c / y);
6644N/A+ c = BASE * (c - y * z.re_fraction[k]);
6644N/A+ }
6644N/A+ if (c < 0)
6644N/A+ {
6644N/A+ mperr ("*** INTEGER OVERFLOW IN DIVIDE_INTEGER, B TOO LARGE ***");
6644N/A+ return new Number.integer (0);
6644N/A+ }
6644N/A+
6644N/A+ mp_normalize (ref z);
6644N/A+ return z;
6644N/A+ }
6644N/A+
6644N/A+ /* HERE NEED SIMULATED DOUBLE-PRECISION DIVISION */
6644N/A+ var j1 = y / BASE;
6644N/A+ var j2 = y - j1 * BASE;
6644N/A+
6644N/A+ /* LOOK FOR FIRST NONZERO DIGIT */
6644N/A+ var c2 = 0;
6644N/A+ do
6644N/A+ {
6644N/A+ c = BASE * c + c2;
6644N/A+ c2 = i < T ? re_fraction[i] : 0;
6644N/A+ i++;
6644N/A+ } while (c < j1 || (c == j1 && c2 < j2));
6644N/A+
6644N/A+ /* COMPUTE T+4 QUOTIENT DIGITS */
6644N/A+ z.re_exponent += (int) (1 - i);
6644N/A+ i--;
6644N/A+
6644N/A+ /* MAIN LOOP FOR LARGE ABS (y) CASE */
6644N/A+ for (var k = 1; k <= T + 4; k++)
6644N/A+ {
6644N/A+ /* GET APPROXIMATE QUOTIENT FIRST */
6644N/A+ var ir = c / (j1 + 1);
6644N/A+
6644N/A+ /* NOW REDUCE SO OVERFLOW DOES NOT OCCUR */
6644N/A+ var iq = c - ir * j1;
6644N/A+ if (iq >= b2)
6644N/A+ {
6644N/A+ /* HERE IQ*B WOULD POSSIBLY OVERFLOW SO INCREASE IR */
6644N/A+ ir++;
6644N/A+ iq -= j1;
6644N/A+ }
6644N/A+
6644N/A+ iq = iq * BASE - ir * j2;
6644N/A+ if (iq < 0)
6644N/A+ {
6644N/A+ /* HERE IQ NEGATIVE SO IR WAS TOO LARGE */
6644N/A+ ir--;
6644N/A+ iq += y;
6644N/A+ }
6644N/A+
6644N/A+ if (i < T)
6644N/A+ iq += re_fraction[i];
6644N/A+ i++;
6644N/A+ var iqj = iq / y;
6644N/A+
6644N/A+ /* r (K) = QUOTIENT, C = REMAINDER */
6644N/A+ z.re_fraction[k - 1] = (int) (iqj + ir);
6644N/A+ c = iq - y * iqj;
6644N/A+
6644N/A+ if (c < 0)
6644N/A+ {
6644N/A+ /* CARRY NEGATIVE SO OVERFLOW MUST HAVE OCCURRED */
6644N/A+ mperr ("*** INTEGER OVERFLOW IN DIVIDE_INTEGER, B TOO LARGE ***");
6644N/A+ return new Number.integer (0);
6644N/A+ }
6644N/A+ }
6644N/A+
6644N/A+ mp_normalize (ref z);
6644N/A+
6644N/A+ /* CARRY NEGATIVE SO OVERFLOW MUST HAVE OCCURRED */
6644N/A+ mperr ("*** INTEGER OVERFLOW IN DIVIDE_INTEGER, B TOO LARGE ***");
6644N/A+ return new Number.integer (0);
6644N/A+ }
6644N/A+
6644N/A+
6644N/A+
6644N/A private Number from_radians (AngleUnit unit)
6644N/A {
6644N/A switch (unit)
6644N/A@@ -1340,23 +2825,146 @@
6644N/A }
6644N/A }
6644N/A
6644N/A+ /* z = sin (x) -1 >= x >= 1, do_sin = 1
6644N/A+ * z = cos (x) -1 >= x >= 1, do_sin = 0
6644N/A+ */
6644N/A+ private Number sin1 (bool do_sin)
6644N/A+ {
6644N/A+ /* sin (0) = 0, cos (0) = 1 */
6644N/A+ if (is_zero ())
6644N/A+ {
6644N/A+ if (do_sin)
6644N/A+ return new Number.integer (0);
6644N/A+ else
6644N/A+ return new Number.integer (1);
6644N/A+ }
6644N/A+
6644N/A+ var t2 = multiply (this);
6644N/A+ if (t2.compare (new Number.integer (1)) > 0)
6644N/A+ mperr ("*** ABS (X) > 1 IN CALL TO SIN1 ***");
6644N/A+
6644N/A+ Number t1;
6644N/A+ int i;
6644N/A+ Number z;
6644N/A+ if (do_sin)
6644N/A+ {
6644N/A+ t1 = this;
6644N/A+ z = t1;
6644N/A+ i = 2;
6644N/A+ }
6644N/A+ else
6644N/A+ {
6644N/A+ t1 = new Number.integer (1);
6644N/A+ z = new Number.integer (0);
6644N/A+ i = 1;
6644N/A+ }
6644N/A+
6644N/A+ /* Taylor series */
6644N/A+ /* POWER SERIES LOOP. REDUCE T IF POSSIBLE */
6644N/A+ var b2 = 2 * int.max (BASE, 64);
6644N/A+ do
6644N/A+ {
6644N/A+ if (T + t1.re_exponent <= 0)
6644N/A+ break;
6644N/A+
6644N/A+ /* IF I*(I+1) IS NOT REPRESENTABLE AS AN INTEGER, THE FOLLOWING
6644N/A+ * DIVISION BY I*(I+1) HAS TO BE SPLIT UP.
6644N/A+ */
6644N/A+ t1 = t2.multiply (t1);
6644N/A+ if (i > b2)
6644N/A+ {
6644N/A+ t1 = t1.divide_integer (-i);
6644N/A+ t1 = t1.divide_integer (i + 1);
6644N/A+ }
6644N/A+ else
6644N/A+ t1 = t1.divide_integer (-i * (i + 1));
6644N/A+ z = t1.add (z);
6644N/A+
6644N/A+ i += 2;
6644N/A+ } while (t1.re_sign != 0);
6644N/A+
6644N/A+ if (!do_sin)
6644N/A+ z = z.add (new Number.integer (1));
6644N/A+
6644N/A+ return z;
6644N/A+ }
6644N/A+
6644N/A+
6644N/A private Number sin_real (AngleUnit unit)
6644N/A {
6644N/A+ /* sin (0) = 0 */
6644N/A+ if (is_zero ())
6644N/A+ return new Number.integer (0);
6644N/A+
6644N/A var x_radians = to_radians (unit);
6644N/A- var z = new Number.integer (0);
6644N/A- var tmp = z.re_num;
6644N/A- tmp.sin (x_radians.re_num, Round.NEAREST);
6644N/A- z.re_num = tmp;
6644N/A+
6644N/A+ var xs = x_radians.re_sign;
6644N/A+ x_radians = x_radians.abs ();
6644N/A+
6644N/A+ /* USE SIN1 IF ABS (X) <= 1 */
6644N/A+ Number z;
6644N/A+ if (x_radians.compare (new Number.integer (1)) <= 0)
6644N/A+ z = x_radians.sin1 (true);
6644N/A+ /* FIND ABS (X) MODULO 2PI */
6644N/A+ else
6644N/A+ {
6644N/A+ z = new Number.pi ().divide_integer (4);
6644N/A+ x_radians = x_radians.divide (z);
6644N/A+ x_radians = x_radians.divide_integer (8);
6644N/A+ x_radians = x_radians.fractional_component ();
6644N/A+
6644N/A+ /* SUBTRACT 1/2, SAVE re_sign AND TAKE ABS */
6644N/A+ x_radians = x_radians.add (new Number.fraction (-1, 2));
6644N/A+ xs = -xs * x_radians.re_sign;
6644N/A+ if (xs == 0)
6644N/A+ return new Number.integer (0);
6644N/A+
6644N/A+ x_radians.re_sign = 1;
6644N/A+ x_radians = x_radians.multiply_integer (4);
6644N/A+
6644N/A+ /* IF NOT LESS THAN 1, SUBTRACT FROM 2 */
6644N/A+ if (x_radians.re_exponent > 0)
6644N/A+ x_radians = x_radians.add (new Number.integer (-2));
6644N/A+
6644N/A+ if (x_radians.is_zero ())
6644N/A+ return new Number.integer (0);
6644N/A+
6644N/A+ x_radians.re_sign = 1;
6644N/A+ x_radians = x_radians.multiply_integer (2);
6644N/A+
6644N/A+ /* NOW REDUCED TO FIRST QUADRANT, IF LESS THAN PI/4 USE
6644N/A+ * POWER SERIES, ELSE COMPUTE COS OF COMPLEMENT
6644N/A+ */
6644N/A+ if (x_radians.re_exponent > 0)
6644N/A+ {
6644N/A+ x_radians = x_radians.add (new Number.integer (-2));
6644N/A+ x_radians = x_radians.multiply (z);
6644N/A+ z = x_radians.sin1 (false);
6644N/A+ }
6644N/A+ else
6644N/A+ {
6644N/A+ x_radians = x_radians.multiply (z);
6644N/A+ z = x_radians.sin1 (true);
6644N/A+ }
6644N/A+ }
6644N/A+
6644N/A+ z.re_sign = xs;
6644N/A return z;
6644N/A }
6644N/A
6644N/A private Number cos_real (AngleUnit unit)
6644N/A {
6644N/A- var x_radians = to_radians (unit);
6644N/A- var tmp = MPFloat.init2 ((Precision) precision);
6644N/A- tmp.cos (x_radians.re_num, Round.NEAREST);
6644N/A- var z = new Number.mpfloat (tmp);
6644N/A- return z;
6644N/A+ /* cos (0) = 1 */
6644N/A+ if (is_zero ())
6644N/A+ return new Number.integer (1);
6644N/A+
6644N/A+ /* Use power series if |x| <= 1 */
6644N/A+ var z = to_radians (unit).abs ();
6644N/A+ if (z.compare (new Number.integer (1)) <= 0)
6644N/A+ return z.sin1 (false);
6644N/A+ else
6644N/A+ /* cos (x) = sin (π/2 - |x|) */
6644N/A+ return new Number.pi ().divide_integer (2).subtract (z).sin (AngleUnit.RADIANS);
6644N/A }
6644N/A
6644N/A private Number bitwise (Number y, BitwiseFunc bitwise_operator, int wordlen)
6644N/A@@ -1370,7 +2978,7 @@
6644N/A offset_out = offset1 > offset2 ? offset1 : offset2;
6644N/A if (offset_out > 0 && (offset_out < offset1 || offset_out < offset2))
6644N/A {
6644N/A- error = ("Overflow. Try a bigger word size");
6644N/A+ mperr ("Overflow. Try a bigger word size");
6644N/A return new Number.integer (0);
6644N/A }
6644N/A
6644N/A@@ -1420,6 +3028,29 @@
6644N/A
6644N/A // FIXME: Re-add overflow and underflow detection
6644N/A
6644N/A+static string? mp_error = null;
6644N/A+
6644N/A+/* THIS ROUTINE IS CALLED WHEN AN ERROR CONDITION IS ENCOUNTERED, AND
6644N/A+ * AFTER A MESSAGE HAS BEEN WRITTEN TO STDERR.
6644N/A+ */
6644N/A+public void mperr (string text)
6644N/A+{
6644N/A+ mp_error = text;
6644N/A+}
6644N/A+
6644N/A+/* Returns error string or null if no error */
6644N/A+// FIXME: Global variable
6644N/A+public string mp_get_error ()
6644N/A+{
6644N/A+ return mp_error;
6644N/A+}
6644N/A+
6644N/A+/* Clear any current error */
6644N/A+public void mp_clear_error ()
6644N/A+{
6644N/A+ mp_error = null;
6644N/A+}
6644N/A+
6644N/A /* Sets z from a string representation in 'text'. */
6644N/A public Number? mp_set_from_string (string str, int default_base = 10)
6644N/A {
6644N/A@@ -1468,7 +3099,6 @@
6644N/A
6644N/A /* Convert integer part */
6644N/A var z = new Number.integer (0);
6644N/A-
6644N/A while (str.get_next_char (ref index, out c))
6644N/A {
6644N/A var i = char_val (c, number_base);
6644N/A@@ -1600,6 +3230,73 @@
6644N/A return null;
6644N/A }
6644N/A
6644N/A+/* RETURNS K = K/GCD AND L = L/GCD, WHERE GCD IS THE
6644N/A+ * GREATEST COMMON DIVISOR OF K AND L.
6644N/A+ * SAVE INPUT PARAMETERS IN LOCAL VARIABLES
6644N/A+ */
6644N/A+public void mp_gcd (ref int64 k, ref int64 l)
6644N/A+{
6644N/A+ var i = k.abs ();
6644N/A+ var j = l.abs ();
6644N/A+ if (j == 0)
6644N/A+ {
6644N/A+ /* IF J = 0 RETURN (1, 0) UNLESS I = 0, THEN (0, 0) */
6644N/A+ k = 1;
6644N/A+ l = 0;
6644N/A+ if (i == 0)
6644N/A+ k = 0;
6644N/A+ return;
6644N/A+ }
6644N/A+
6644N/A+ /* EUCLIDEAN ALGORITHM LOOP */
6644N/A+ do
6644N/A+ {
6644N/A+ i %= j;
6644N/A+ if (i == 0)
6644N/A+ {
6644N/A+ k = k / j;
6644N/A+ l = l / j;
6644N/A+ return;
6644N/A+ }
6644N/A+ j %= i;
6644N/A+ } while (j != 0);
6644N/A+
6644N/A+ /* HERE J IS THE GCD OF K AND L */
6644N/A+ k = k / i;
6644N/A+ l = l / i;
6644N/A+}
6644N/A+
6644N/A+// FIXME: Is r.re_fraction large enough? It seems to be in practise but it may be T+4 instead of T
6644N/A+// FIXME: There is some sort of stack corruption/use of unitialised variables here. Some functions are
6644N/A+// using stack variables as x otherwise there are corruption errors. e.g. "Cos (45) - 1/Sqrt (2) = -0"
6644N/A+// (try in scientific mode)
6644N/A+public void mp_normalize (ref Number x)
6644N/A+{
6644N/A+ int start_index;
6644N/A+
6644N/A+ /* Find first non-zero digit */
6644N/A+ for (start_index = 0; start_index < SIZE && x.re_fraction[start_index] == 0; start_index++);
6644N/A+
6644N/A+ /* Mark as zero */
6644N/A+ if (start_index >= SIZE)
6644N/A+ {
6644N/A+ x.re_sign = 0;
6644N/A+ x.re_exponent = 0;
6644N/A+ return;
6644N/A+ }
6644N/A+
6644N/A+ /* Shift left so first digit is non-zero */
6644N/A+ if (start_index > 0)
6644N/A+ {
6644N/A+ x.re_exponent -= start_index;
6644N/A+ var i = 0;
6644N/A+ for (; (i + start_index) < SIZE; i++)
6644N/A+ x.re_fraction[i] = x.re_fraction[i + start_index];
6644N/A+ for (; i < SIZE; i++)
6644N/A+ x.re_fraction[i] = 0;
6644N/A+ }
6644N/A+}
6644N/A+
6644N/A /* Returns true if x is cannot be represented in a binary word of length 'wordlen' */
6644N/A public bool mp_is_overflow (Number x, int wordlen)
6644N/A {
6644N/A
6644N/A--- gnome-calculator-3.18.2/configure.ac 2016-08-09 14:02:40.743464258 +0000
6644N/A+++ gnome-calculator-3.18.2/configure.ac 2016-08-09 14:02:52.623794159 +0000
6644N/A@@ -25,9 +25,6 @@
6644N/A
6644N/A AC_SUBST([GLIB_REQUIRED])
6644N/A
6644N/A-# FIXME We link to it too, so this check needs to be better....
6644N/A-AC_CHECK_HEADER([mpfr.h], [], [AC_MSG_ERROR([The mpfr header is missing. Please, install mpfr])])
6644N/A-
6644N/A PKG_CHECK_MODULES(LIBCALCULATOR, [
6644N/A glib-2.0 >= $GLIB_REQUIRED
6644N/A gio-2.0 >= $GLIB_REQUIRED