Double Precision Real Numbers

EXAMPLES:

We create the real double vector space of dimension 3:

sage: V = RDF^3; V
Vector space of dimension 3 over Real Double Field

Notice that this space is unique.

sage: V is RDF^3
True
sage: V is FreeModule(RDF, 3)
True
sage: V is VectorSpace(RDF, 3)
True

Also, you can instantly create a space of large dimension.

sage: V = RDF^10000

TESTS:

Test NumPy conversions::

    sage: RDF(1).__array_interface__
    {'typestr': '=f8'}
    sage: import numpy
    sage: numpy.array([RDF.pi()]).dtype
    dtype('float64')
class sage.rings.real_double.RealDoubleElement

An approximation to a real number using double precision floating point numbers. Answers derived from calculations with such approximations may differ from what they would be if those calculations were performed with true real numbers. This is due to the rounding errors inherent to finite precision calculations.

NaN()

EXAMPLES:

sage: RDF.NaN()
NaN
__abs__()

Returns the absolute value of self.

EXAMPLES:

sage: abs(RDF(1.5))
1.5
sage: abs(RDF(-1.5))
1.5
__complex__()

EXAMPLES:

sage: a = 2303
sage: RDF(a)
2303.0
sage: complex(RDF(a))
(2303+0j)
__copy__()

Return copy of self, which since self is immutable, is just self.

EXAMPLES:

sage: r = RDF('-1.6')
sage: r.__copy__() is r
True
__eq__()
x.__eq__(y) <==> x==y
__float__()

Return self as a python float.

EXAMPLES:

sage: float(RDF(1.5))
1.5
sage: type(float(RDF(1.5)))
<type 'float'>
__ge__()
x.__ge__(y) <==> x>=y
__gt__()
x.__gt__(y) <==> x>y
__hash__()
x.__hash__() <==> hash(x)
__init__()
x.__init__(...) initializes x; see x.__class__.__doc__ for signature
__int__()

Returns integer truncation of this real number.

EXAMPLES:

sage: int(RDF(2.99))
2
sage: int(RDF(-2.99))
-2
__invert__()

Compute the multiplicative inverse of self.

EXAMPLES:

sage: a = RDF(-1.5)*RDF(2.5)
sage: a.__invert__()
-0.266666666667
sage: ~a
-0.266666666667
__le__()
x.__le__(y) <==> x<=y
__long__()

Returns long integer truncation of this real number.

EXAMPLES:

sage: int(RDF(10e15))
10000000000000000L                   # 32-bit
10000000000000000                    # 64-bit
sage: long(RDF(2^100)) == 2^100
True
__lshift__()
LShifting a double is not supported; nor is lshifting a RealDoubleElement.
__lt__()
x.__lt__(y) <==> x<y
__ne__()
x.__ne__(y) <==> x!=y
__neg__()

Negates a real number.

EXAMPLES:

sage: -RDF('-1.5')
1.5
static __new__()
T.__new__(S, ...) -> a new object with type S, a subtype of T
__pow__()
x.__pow__(y[, z]) <==> pow(x, y[, z])
__reduce__()

EXAMPLES:

sage: a = RDF(-2.7)
sage: loads(dumps(a)) == a
True
__repr__()

Return print version of self.

EXAMPLES:

sage: a = RDF(2); a
2.0
sage: a^2
4.0
__rlshift__()
x.__rlshift__(y) <==> y<<x
__rpow__()
y.__rpow__(x[, z]) <==> pow(x, y[, z])
__rrshift__()
x.__rrshift__(y) <==> y>>x
__rshift__()
RShifting a double is not supported; nor is rshifting a RealDoubleElement.
_add_()

Add two real numbers with the same parent.

EXAMPLES:

sage: RDF('-1.5') + RDF('2.5')
1.0
_complex_double_()

EXAMPLES:

sage: CDF(RDF(1/3))
0.333333333333
_complex_mpfr_field_()

EXAMPLES:

sage: a = RDF(1/3)
sage: CC(a)
0.333333333333333
sage: a._complex_mpfr_field_(CC)
0.333333333333333

If we coerce to a higher-precision field the extra bits appear random; they are actually 0’s in base 2.

sage: a._complex_mpfr_field_(ComplexField(100))
0.33333333333333331482961625625
sage: a._complex_mpfr_field_(ComplexField(100)).str(2)
'0.01010101010101010101010101010101010101010101010101010100000000000000000000000000000000000000000000000'
_div_()

EXAMPLES:

sage: RDF('-1.5') / RDF('2.5')
-0.6
sage: RDF(1)/RDF(0)
+infinity
_iadd_()

EXAMPLES:

sage: a = RDF(0.5)
sage: a += RDF(3); a
3.5
_idiv_()

EXAMPLES:

sage: a = RDF(1.5)
sage: a /= RDF(2); a
0.75
sage: a /= RDF(0); a
+infinity
_im_gens_()
_imul_()

EXAMPLES:

sage: a = RDF(2.5)
sage: a *= RDF(3); a
7.5
_integer_()

If this floating-point number is actually an integer, return that integer. Otherwise, raise an exception.

EXAMPLES:

sage: ZZ(RDF(237.0))
237
sage: ZZ(RDF(0.0/0.0))
...
ValueError: cannot convert float NaN to integer
sage: ZZ(RDF(1.0/0.0))
...
OverflowError: cannot convert float infinity to integer
sage: ZZ(RDF(-123456789.0))
-123456789
sage: ZZ(RDF((2.0))^290)
1989292945639146568621528992587283360401824603189390869761855907572637988050133502132224
sage: ZZ(RDF(-2345.67))
...
TypeError: Cannot convert non-integral float to integer
_interface_init_()

Returns self formatted as a string, suitable as input to another computer algebra system. (This is the default function used for exporting to other computer algebra systems.)

EXAMPLES:

sage: s1 = RDF(sin(1)); s1
0.841470984808
sage: s1._interface_init_()
'0.8414709848078965'
sage: s1 == RDF(gp(s1))
True
_isub_()

EXAMPLES:

sage: a = RDF(0.5)
sage: a -= RDF(3); a
-2.5
_latex_()

EXAMPLE:

sage: RDF(2e-100)._latex_()
'2 \times 10^{-100}'
_magma_init_()

Return a string representation of self in the Magma language.

EXAMPLES:

sage: RDF(10.5)
10.5
sage: magma(RDF(10.5)) # optional - magma
10.5000000000000
_mul_()

Multiply two real numbers with the same parent.

EXAMPLES:

sage: RDF('-1.5') * RDF('2.5')
-3.75
_pari_()

EXAMPLES:

sage: RDF(1.5)._pari_()
1.50000000000000
_rpy_()

Returns self.__float__() for rpy to convert into the appropriate R object.

EXAMPLES:

sage: n = RDF(2.0)
sage: n._rpy_()
2.0
sage: type(n._rpy_())
<type 'float'>
_sage_input_()

Produce an expression which will reproduce this value when evaluated.

EXAMPLES:

sage: sage_input(RDF(NaN), verify=True)
# Verified
RDF(NaN)
sage: sage_input(RDF(-infinity), verify=True)
# Verified
-RDF(infinity)
sage: sage_input(RDF(-infinity)*polygen(RDF), verify=True)
# Verified
R.<x> = RDF[]
-RDF(infinity)*x + RDF(NaN)
sage: sage_input(RDF(pi), verify=True)
# Verified
RDF(3.1415926535897931)
sage: sage_input(RDF(-e), verify=True, preparse=False)
# Verified
-RDF(2.718281828459045...)
sage: sage_input(RDF(pi)*polygen(RDF), verify=True, preparse=None)
# Verified
R = RDF['x']
x = R.gen()
3.1415926535897931*x
sage: from sage.misc.sage_input import SageInputBuilder
sage: sib = SageInputBuilder()
sage: RDF(22/7)._sage_input_(sib, True)
{atomic:3.1428571428571428}
sage: RDF(22/7)._sage_input_(sib, False)
{call: {atomic:RDF}({atomic:3.1428571428571428})}
_sub_()

Subtract two real numbers with the same parent.

EXAMPLES:

sage: RDF('-1.5') - RDF('2.5')
-4.0
abs()

Returns the absolute value of self.

EXAMPLES:

sage: RDF(1e10).abs()
10000000000.0
sage: RDF(-1e10).abs()
10000000000.0
acosh()

Returns the hyperbolic inverse cosine of this number

EXAMPLES:

sage: q = RDF.pi()/2 
sage: i = q.cosh() ; i
2.50917847866
sage: abs(i.acosh()-q) < 1e-15
True
agm()

Return the arithmetic-geometric mean of self and other. The arithmetic-geometric mean is the common limit of the sequences u_n and v_n, where u_0 is self, v_0 is other, u_{n+1} is the arithmetic mean of u_n and v_n, and v_{n+1} is the geometric mean of u_n and v_n. If any operand is negative, the return value is NaN.

EXAMPLES:

sage: a = RDF(1.5)
sage: b = RDF(2.3)
sage: a.agm(b)
1.87864845581

The arithmetic-geometric mean always lies between the geometric and arithmetic mean.

sage: sqrt(a*b) < a.agm(b) < (a+b)/2
True
algdep()

Returns a polynomial of degree at most n which is approximately satisfied by this number. Note that the returned polynomial need not be irreducible, and indeed usually won’t be if this number is a good approximation to an algebraic number of degree less than n.

ALGORITHM: Uses the PARI C-library algdep command.

EXAMPLE:

sage: r = RDF(2).sqrt(); r
1.41421356237
sage: r.algdep(5)
x^2 - 2
algebraic_dependency()

Returns a polynomial of degree at most n which is approximately satisfied by this number. Note that the returned polynomial need not be irreducible, and indeed usually won’t be if this number is a good approximation to an algebraic number of degree less than n.

ALGORITHM: Uses the PARI C-library algdep command.

EXAMPLE:

sage: r = sqrt(RDF(2)); r
1.41421356237
sage: r.algdep(5)
x^2 - 2
arccos()

Returns the inverse cosine of this number

EXAMPLES:

sage: q = RDF.pi()/3
sage: i = q.cos()
sage: i.arccos() == q
True
arcsin()

Returns the inverse sine of this number

EXAMPLES:

sage: q = RDF.pi()/5
sage: i = q.sin()
sage: i.arcsin() == q
True
arcsinh()

Returns the hyperbolic inverse sine of this number

EXAMPLES:

sage: q = RDF.pi()/2 
sage: i = q.sinh() ; i
2.30129890231
sage: abs(i.arcsinh()-q) < 1e-15
True
arctan()

Returns the inverse tangent of this number

EXAMPLES:

sage: q = RDF.pi()/5
sage: i = q.tan()
sage: i.arctan() == q
True
arctanh()

Returns the hyperbolic inverse tangent of this number

EXAMPLES:

sage: q = RDF.pi()/2 
sage: i = q.tanh() ; i
0.917152335667
sage: i.arctanh() - q      # output is random, depending on arch.
-4.4408920985e-16
ceil()

Returns the ceiling of this number

OUTPUT: integer

EXAMPLES:

sage: RDF(2.99).ceil()
3
sage: RDF(2.00).ceil()
2
sage: RDF(-5/2).ceil()
-2
ceiling()

Returns the ceiling of this number

OUTPUT: integer

EXAMPLES:

sage: RDF(2.99).ceil()
3
sage: RDF(2.00).ceil()
2
sage: RDF(-5/2).ceil()
-2
cos()

Returns the cosine of this number

EXAMPLES:

sage: t=RDF.pi()/2
sage: t.cos()
6.12323399574e-17
cosh()

Returns the hyperbolic cosine of this number

EXAMPLES:

sage: q = RDF.pi()/12
sage: q.cosh()
1.0344656401
coth()

This function returns the hyperbolic cotangent.

EXAMPLES:

sage: RDF(pi).coth()
1.0037418732
sage: CDF(pi).coth()
1.0037418732
csch()

This function returns the hyperbolic cosecant.

EXAMPLES:

sage: RDF(pi).csch()
0.08658953753
sage: CDF(pi).csch()
0.08658953753
cube_root()

Return the cubic root (defined over the real numbers) of self.

EXAMPLES:

sage: r = RDF(125.0); r.cube_root()
5.0
sage: r = RDF(-119.0)
sage: r.cube_root()^3 - r         # output is random, depending on arch. 
0.0
erf()

Returns the value of the error function on self.

EXAMPLES:

sage: RDF(6).erf()
1.0
exp()

Returns e^\mathtt{self}

EXAMPLES:

sage: r = RDF(0.0)
sage: r.exp()
1.0
sage: r = RDF('32.3')
sage: a = r.exp(); a
1.06588847275e+14
sage: a.log()
32.3
sage: r = RDF('-32.3')
sage: r.exp()
9.3818445885e-15
sage: RDF(1000).exp()
+infinity
exp10()

Returns 10^\mathtt{self}

EXAMPLES:

sage: r = RDF(0.0)
sage: r.exp10()
1.0
sage: r = RDF(32.0)
sage: r.exp10()
1e+32
sage: r = RDF(-32.3)
sage: r.exp10()
5.01187233627e-33
exp2()

Returns 2^\mathtt{self}

EXAMPLES:

sage: r = RDF(0.0)
sage: r.exp2()
1.0
sage: r = RDF(32.0)
sage: r.exp2()
4294967296.0
sage: r = RDF(-32.3)
sage: r.exp2()
1.89117248253e-10
floor()

Returns the floor of this number

EXAMPLES:

sage: RDF(2.99).floor()
2
sage: RDF(2.00).floor()
2
sage: RDF(-5/2).floor()
-3
frac()

frac returns a real number > -1 and < 1. it satisfies the relation: x = x.trunc() + x.frac()

EXAMPLES:

sage: RDF(2.99).frac()
0.99
sage: RDF(2.50).frac()
0.5
sage: RDF(-2.79).frac()
-0.79
gamma()

The Euler gamma function. Return gamma of self.

EXAMPLES:

sage: RDF(6).gamma()
120.0
sage: RDF(1.5).gamma()
0.886226925453
hypot()

Computes the value \sqrt(self^2 + other^2) in such a way as to avoid overflow.

EXAMPLES:

sage: x = RDF(4e300); y = RDF(3e300); 
sage: x.hypot(y)
5e+300
sage: sqrt(x^2+y^2) # overflow
+infinity
imag()

Returns the imaginary part of this number. (hint: it’s zero.)

EXAMPLES:

sage: a = RDF(3)
sage: a.imag()
0.0
integer_part()

If in decimal this number is written n.defg, returns n.

EXAMPLES:

sage: r = RDF('-1.6')
sage: a = r.integer_part(); a
-1
sage: type(a)
<type 'sage.rings.integer.Integer'>
sage: r = RDF(0.0/0.0)
sage: a = r.integer_part()
...
TypeError: Attempt to get integer part of NaN
is_NaN()

EXAMPLES:

sage: RDF(1).is_NaN()
False
sage: a = RDF(0)/RDF(0)
sage: a.is_NaN()
True
is_infinity()

EXAMPLES:

sage: a = RDF(2); b = RDF(0)
sage: (a/b).is_infinity()
True
sage: (b/a).is_infinity()
False
is_negative_infinity()

EXAMPLES:

sage: a = RDF(2)/RDF(0)
sage: a.is_negative_infinity()
False
sage: a = RDF(-3)/RDF(0)
sage: a.is_negative_infinity()
True
is_positive_infinity()

EXAMPLES:

sage: a = RDF(1)/RDF(0)
sage: a.is_positive_infinity()
True
sage: a = RDF(-1)/RDF(0)
sage: a.is_positive_infinity()
False
is_square()

Returns whether or not this number is a square in this field. For the real numbers, this is True if and only if self is non-negative.

EXAMPLES:

sage: RDF(3.5).is_square()
True
sage: RDF(0).is_square()
True
sage: RDF(-4).is_square()
False
log()

EXAMPLES:

sage: RDF(2).log()
0.69314718056
sage: RDF(2).log(2)
1.0
sage: RDF(2).log(pi)
0.605511561398
sage: RDF(2).log(10)
0.301029995664
sage: RDF(2).log(1.5)
1.70951129135
sage: RDF(0).log()
-infinity
sage: RDF(-1).log()
NaN
log10()

Returns log to the base 10 of self

EXAMPLES:

sage: r = RDF('16.0'); r.log10()
1.20411998266
sage: r.log() / RDF(log(10))
1.20411998266
sage: r = RDF('39.9'); r.log10()
1.60097289569
log2()

Returns log to the base 2 of self

EXAMPLES:

sage: r = RDF(16.0)
sage: r.log2()
4.0
sage: r = RDF(31.9); r.log2()
4.99548451888
logpi()

Returns log to the base pi of self

EXAMPLES:

sage: r = RDF(16); r.logpi()
2.42204624559
sage: r.log() / RDF(log(pi))
2.42204624559
sage: r = RDF('39.9'); r.logpi()
3.22030233461
multiplicative_order()

Returns n such that self^n == 1.

Only \pm 1 have finite multiplicative order.

EXAMPLES:

sage: RDF(1).multiplicative_order()
1
sage: RDF(-1).multiplicative_order()
2
sage: RDF(3).multiplicative_order()
+Infinity
nan()

EXAMPLES:

sage: RDF.nan()
NaN
nth_root()

Returns the n^{th} root of self.

INPUT:

  • n - an integer

OUTPUT: an real or complex double

The output is complex if self is negative and n is even.

EXAMPLES:

sage: r = RDF(-125.0); r.nth_root(3)
-5.0
sage: r.nth_root(5)
-2.6265278044
sage: RDF(-2).nth_root(5)^5
-2.0
sage: RDF(-1).nth_root(5)^5
-1.0
sage: RDF(3).nth_root(10)^10
3.0
sage: RDF(-1).nth_root(2)
6.12323399574e-17 + 1.0*I
sage: RDF(-1).nth_root(4)
0.707106781187 + 0.707106781187*I
parent()

Return the real double field, which is the parent of self.

EXAMPLES:

sage: a = RDF(2.3)
sage: a.parent()
Real Double Field
sage: parent(a)
Real Double Field
prec()

Returns the precision of this number in bits.

Always returns 53.

EXAMPLES:

sage: RDF(0).prec()
53
real()

Returns itself - we’re already real.

EXAMPLES:

sage: a = RDF(3)
sage: a.real()
3.0
restrict_angle()

Returns a number congruent to self mod 2\pi that lies in the interval (-\pi, \pi].

Specifically, it is the unique x \in (-\pi, \pi] such that `self = x + 2pi n` for some n \in \ZZ.

EXAMPLES:

sage: RDF(pi).restrict_angle()
3.14159265359
sage: RDF(pi + 1e-10).restrict_angle()
-3.14159265349
sage: RDF(1+10^10*pi).restrict_angle()
0.9999977606...
round()

Given real number x, rounds up if fractional part is greater than .5, rounds down if fractional part is less than .5.

EXAMPLES:

sage: RDF(0.49).round()
0
sage: a=RDF(0.51).round(); a
1
sech()

This function returns the hyperbolic secant.

EXAMPLES:

sage: RDF(pi).sech()
0.0862667383341
sage: CDF(pi).sech()
0.0862667383341
sign()

Returns -1,0, or 1 if self is negative, zero, or positive; respectively.

EXAMPLES:

sage: RDF(-1.5).sign()
-1
sage: RDF(0).sign()
0
sage: RDF(2.5).sign()
1
sin()

Returns the sine of this number

EXAMPLES:

sage: RDF(2).sin()
0.909297426826
sincos()

Returns a pair consisting of the sine and cosine.

EXAMPLES:

sage: t = RDF.pi()/6
sage: t.sincos()
(0.5, 0.866025403784)
sinh()

Returns the hyperbolic sine of this number

EXAMPLES:

sage: q = RDF.pi()/12
sage: q.sinh()
0.264800227602
sqrt()

The square root function.

INPUT:

  • extend - bool (default: True); if True, return a square root in a complex field if necessary if self is negative; otherwise raise a ValueError
  • all - bool (default: False); if True, return a list of all square roots.

EXAMPLES:

sage: r = RDF(4.0)
sage: r.sqrt()
2.0
sage: r.sqrt()^2 == r
True
sage: r = RDF(4344)
sage: r.sqrt()
65.9090282131
sage: r.sqrt()^2 - r
0.0
sage: r = RDF(-2.0)
sage: r.sqrt()
1.41421356237*I
sage: RDF(2).sqrt(all=True)
[1.41421356237, -1.41421356237]
sage: RDF(0).sqrt(all=True)
[0.0]
sage: RDF(-2).sqrt(all=True)
[1.41421356237*I, -1.41421356237*I]
str()

Return string representation of self.

EXAMPLES:

sage: a = RDF('4.5'); a.str()
'4.5'
sage: a = RDF('49203480923840.2923904823048'); a.str()
'4.92034809238e+13'
sage: a = RDF(1)/RDF(0); a.str()
'+infinity'
sage: a = -RDF(1)/RDF(0); a.str()
'-infinity'
sage: a = RDF(0)/RDF(0); a.str()
'NaN'

We verify consistency with RR (mpfr reals):

sage: str(RR(RDF(1)/RDF(0))) == str(RDF(1)/RDF(0))
True
sage: str(RR(-RDF(1)/RDF(0))) == str(-RDF(1)/RDF(0))
True
sage: str(RR(RDF(0)/RDF(0))) == str(RDF(0)/RDF(0))
True
tan()

Returns the tangent of this number

EXAMPLES:

sage: q = RDF.pi()/3
sage: q.tan()
1.73205080757
sage: q = RDF.pi()/6
sage: q.tan()
0.57735026919
tanh()

Returns the hyperbolic tangent of this number

EXAMPLES:

sage: q = RDF.pi()/12
sage: q.tanh()
0.255977789246
trunc()

Truncates this number (returns integer part).

EXAMPLES:

sage: RDF(2.99).trunc()
2
sage: RDF(-2.00).trunc()
-2
sage: RDF(0.00).trunc()
0
ulp()

Returns the unit of least precision of self, which is the weight of the least significant bit of self. Unless self is exactly a power of two, it is gap between this number and the next closest distinct number that can be represented.

EXAMPLES:

sage: a = RDF(1)
sage: a - a.ulp() == a
False
sage: a - a.ulp()/2 == a
True

sage: a = RDF.pi()
sage: b = a + a.ulp()
sage: (a+b)/2 in [a,b]
True

sage: a = RDF(1)/RDF(0); a
+infinity
sage: a.ulp()
+infinity
sage: (-a).ulp()
+infinity
sage: a = RR('nan')
sage: a.ulp() is a
True
zeta()

Return the Riemann zeta function evaluated at this real number.

Note

PARI is vastly more efficient at computing the Riemann zeta function. See the example below for how to use it.

EXAMPLES:

sage: RDF(2).zeta()
1.64493406685
sage: RDF.pi()^2/6
1.64493406685
sage: RDF(-2).zeta()       # slightly random-ish arch dependent output
-2.37378795339e-18
sage: RDF(1).zeta()
+infinity
sage.rings.real_double.RealDoubleField()

Return the unique instance of the Real Double Field.

EXAMPLES:

sage: RealDoubleField() is RealDoubleField()
True
class sage.rings.real_double.RealDoubleField_class

An approximation to the field of real numbers using double precision floating point numbers. Answers derived from calculations in this approximation may differ from what they would be if those calculations were performed in the true field of real numbers. This is due to the rounding errors inherent to finite precision calculations.

EXAMPLES:

sage: RR == RDF
False
sage: RDF == RealDoubleField()    # RDF is the shorthand
True
sage: RDF(1)
1.0
sage: RDF(2/3)
0.666666666667

A TypeError is raised if the coercion doesn’t make sense:

sage: RDF(QQ['x'].0)
...
TypeError: cannot coerce nonconstant polynomial to float
sage: RDF(QQ['x'](3))
3.0

One can convert back and forth between double precision real numbers and higher-precision ones, though of course there may be loss of precision:

sage: a = RealField(200)(2).sqrt(); a
1.4142135623730950488016887242096980785696718753769480731767
sage: b = RDF(a); b
1.41421356237
sage: a.parent()(b)
1.4142135623730951454746218587388284504413604736328125000000
sage: a.parent()(b) == b
True
sage: b == RR(a)
True
NaN()

EXAMPLES:

sage: RDF.NaN()
NaN
__cmp__()
x.__cmp__(y) <==> cmp(x,y)
__hash__()
x.__hash__() <==> hash(x)
__init__()
x.__init__(...) initializes x; see x.__class__.__doc__ for signature
static __new__()
T.__new__(S, ...) -> a new object with type S, a subtype of T
__reduce__()

EXAMPLES:

sage: loads(dumps(RDF)) is RDF
True
__repr__()

Print out this real double field.

EXAMPLES:

sage: RealDoubleField()
Real Double Field
sage: RDF
Real Double Field
_coerce_map_from_()

Canonical coercion of x to the real double field.

The rings that canonically coerce to the real double field are:

  • the real double field itself
  • int, long, integer, and rational rings
  • real mathematical constants
  • the MPFR real field with <= 53 bits of precision

EXAMPLES:

sage: RDF.coerce(5)
5.0
sage: RDF.coerce(9499294r)
9499294.0
sage: RDF.coerce(61/3)
20.3333333333
sage: parent(RDF(3) + CDF(5))
Complex Double Field
sage: parent(CDF(5) + RDF(3))
Complex Double Field
sage: CDF.gen(0) + 5.0
5.0 + 1.0*I
sage: RLF(2/3) + RDF(1)
1.66666666667
_latex_()
_magma_init_()

Return a string representation of self in the Magma language.

EXAMPLES:

Magma handles precision in decimal digits, so we lose a bit:

sage: magma(RDF) # optional - magma
Real field of precision 15
sage: 10^15 < 2^53 < 10^16
True

When we convert back from Magma, we convert to a generic real field that has 53 bits of precision:

sage: magma(RDF).sage() # optional - magma
Real Field with 53 bits of precision
_sage_input_()

Produce an expression which will reproduce this value when evaluated.

EXAMPLES:
sage: sage_input(RDF, verify=True) # Verified RDF sage: from sage.misc.sage_input import SageInputBuilder sage: RDF._sage_input_(SageInputBuilder(), False) {atomic:RDF}
algebraic_closure()

Returns the algebraic closure of self, ie, the complex double field.

EXAMPLES:

sage: RDF.algebraic_closure()
Complex Double Field
characteristic()

Returns 0, since the field of real numbers has characteristic 0.

EXAMPLES:

sage: RDF.characteristic()
0
complex_field()

Returns the complex field with the same precision as self, ie, the complex double field.

EXAMPLES:

sage: RDF.complex_field()
Complex Double Field
construction()

Returns the functorial construction of self, namely, completion of the rational numbers with respect to the prime at \infty.

Also preserves other information that makes this field unique (i.e. the Real Double Field).

EXAMPLES:

sage: c, S = RDF.construction(); S
Rational Field
sage: RDF == c(S)
True
euler_constant()

Returns Euler’s gamma constant to double precision

EXAMPLES:

sage: RDF.euler_constant()
0.577215664902
factorial()

Return the factorial of the integer n as a real number.

EXAMPLES:

sage: RDF.factorial(100)
9.33262154439e+157
gen()

Return the generator of the real double field.

EXAMPLES:

sage: RDF.0
1.0
sage: RDF.gens()
(1.0,)
is_atomic_repr()

Returns True, to signify that elements of this field print without sums, so parenthesis aren’t required, e.g., in coefficients of polynomials.

EXAMPLES:

sage: RDF.is_atomic_repr()
True
is_exact()

Returns False, because doubles are not exact.

EXAMPLE:

sage: RDF.is_exact()
False
is_finite()

Returns False, since the field of real numbers is not finite. Technical note: There exists an upper bound on the double representation.

EXAMPLES:

sage: RDF.is_finite()
False
log2()

Returns log(2) to the precision of this field.

EXAMPLES:

sage: RDF.log2() 
0.69314718056
sage: RDF(2).log()
0.69314718056
name()
nan()

EXAMPLES:

sage: RDF.nan()
NaN
ngens()
pi()

Returns pi to double-precision.

EXAMPLES:

sage: RDF.pi()
3.14159265359
sage: RDF.pi().sqrt()/2
0.886226925453
prec()

Return the precision of this real double field in bits.

Always returns 53.

EXAMPLES:

sage: RDF.prec()
53
random_element()

Return a random element of this real double field in the interval [min, max].

EXAMPLES:

sage: RDF.random_element()
0.736945423566
sage: RDF.random_element(min=100, max=110)
102.815947352
to_prec()

Returns the real field to the specified precision. As doubles have fixed precision, this will only return a real double field if prec is exactly 53.

EXAMPLES:

sage: RDF.to_prec(52)
Real Field with 52 bits of precision
sage: RDF.to_prec(53)
Real Double Field
zeta()

Return an n-th root of unity in the real field, if one exists, or raise a ValueError otherwise.

EXAMPLES:

sage: RDF.zeta()
-1.0
sage: RDF.zeta(1)
1.0
sage: RDF.zeta(5)
...
ValueError: No 5th root of unity in self
class sage.rings.real_double.ToRDF
__init__()
x.__init__(...) initializes x; see x.__class__.__doc__ for signature
static __new__()
T.__new__(S, ...) -> a new object with type S, a subtype of T
_call_()
_repr_type()
sage.rings.real_double.is_RealDoubleElement()
sage.rings.real_double.is_RealDoubleField()

Returns True if x is the field of real double precision numbers.

EXAMPLE:

sage: from sage.rings.real_double import is_RealDoubleField
sage: is_RealDoubleField(RDF)
True
sage: is_RealDoubleField(RealField(53))
False
sage.rings.real_double.pool_stats()
sage.rings.real_double.time_alloc()
sage.rings.real_double.time_alloc_list()

Previous topic

Fixed and Arbitrary Precision Numerical Fields

Next topic

Double Precision Complex Numbers

This Page