Opened 5 years ago

Closed 5 years ago

#17285 closed defect (fixed)

CIF is missing many functions

Reported by: jdemeyer Owned by:
Priority: major Milestone: sage-6.6
Component: basic arithmetic Keywords:
Cc: tmonteil, vdelecroix Merged in:
Authors: Vincent Delecroix Reviewers: Jeroen Demeyer
Report Upstream: N/A Work issues:
Branch: afccfe6 (Commits) Commit: afccfe62eac92dd2cea900f3bac44600e633aa7d
Dependencies: #17130 Stopgaps:

Description (last modified by vdelecroix)

sage: CIF(cos(3/2))
---------------------------------------------------------------------------
TypeError                                 Traceback (most recent call last)
<ipython-input-31-4eae9038f0b3> in <module>()
----> 1 CIF(cos(Integer(3)/Integer(2)))

/usr/local/src/sage-config/local/lib/python2.7/site-packages/sage/rings/complex_interval_field.pyc in __call__(self, x, im)
    378 
    379             try:
--> 380                 return x._complex_mpfi_( self )
    381             except AttributeError:
    382                 pass

/usr/local/src/sage-config/local/lib/python2.7/site-packages/sage/symbolic/expression.so in sage.symbolic.expression.Expression._complex_mpfi_ (build/cythonized/sage/symbolic/expression.cpp:8043)()

The reason that this fails is:

sage: CIF(3/2).cos()
---------------------------------------------------------------------------
AttributeError                            Traceback (most recent call last)
<ipython-input-13-96dae572cb36> in <module>()
----> 1 CIF(Integer(3)/Integer(2)).cos()

/usr/local/src/sage-config/local/lib/python2.7/site-packages/sage/structure/element.so in sage.structure.element.Element.__getattr__ (build/cythonized/sage/structure/element.c:4068)()

/usr/local/src/sage-config/local/lib/python2.7/site-packages/sage/structure/misc.so in sage.structure.misc.getattr_from_other_class (build/cythonized/sage/structure/misc.c:1631)()

AttributeError: 'sage.rings.complex_interval.ComplexIntervalFieldElement' object has no attribute 'cos'

This problem was also encoutered on this ask thread.

follow up: #18135 and #18136

Change History (22)

comment:1 Changed 5 years ago by jdemeyer

  • Component changed from symbolics to basic arithmetic
  • Description modified (diff)
  • Summary changed from CIF(cos(3/2)) fails to CIF is missing many functions

comment:2 follow-up: Changed 5 years ago by tmonteil

  • Cc tmonteil added
  • Description modified (diff)

Are you working on this ?

Is a solution involving the use of the (already defined) .exp() method acceptable on the short term (it may not lead to the thinest possible intervals), or do you think about some more accurate way of writing those ?

comment:3 in reply to: ↑ 2 Changed 5 years ago by jdemeyer

Replying to tmonteil:

Are you working on this ?

No.

Is a solution involving the use of the (already defined) .exp() method acceptable on the short term (it may not lead to the thinest possible intervals), or do you think about some more accurate way of writing those ?

I think it's acceptable to rely on exp() but instead, I would rather use the structure of the exp() method: completely reduce the computation to a computation in RIF.

comment:4 Changed 5 years ago by fredrik.johansson

For best accuracy, one should use cos(a+bi) = cos(a)*cosh(b) + sin(a)*sinh(b)*i and sin(a+bi) = sin(a)*cosh(b) + cos(a)*sinh(b)*i (if I got those formulas right).

There's a mpfi_sin_cos but it looks like there's no mpfi_sinh_cosh yet.

comment:5 Changed 5 years ago by vdelecroix

  • Cc vdelecroix added

This is also responsible for

sage: a = RLF(cos(2/3))
sage: b = CLF(QQbar(-2).sqrt())
sage: a + b
<repr(<sage.rings.real_lazy.LazyBinop at 0x7fb7ffc28680>) failed: TypeError: unable to simplify to a complex interval approximation>

comment:6 Changed 5 years ago by vdelecroix

  • Authors set to Vincent Delecroix
  • Branch set to u/vdelecroix/17285
  • Commit set to dd4cf4e61cbc44485985e1f74de201c87ad06bc9
  • Status changed from new to needs_review

New commits:

dd4cf4eTrac 17285: CIF trigonometric functions

comment:7 Changed 5 years ago by vdelecroix

This is not the fastest and safest way possible (the obtained interval is sometimes too large). But at least, the given implementation has the feature that trigonometric functions preserve real and imaginary axes...

Vincent

comment:8 Changed 5 years ago by vdelecroix

Note: the link to the ask question is dead!

comment:9 follow-up: Changed 5 years ago by jdemeyer

I don't like the added complexity of having 3 branches in every function. Why not just use the generic case everywhere?

comment:10 in reply to: ↑ 9 Changed 5 years ago by vdelecroix

Replying to jdemeyer:

I don't like the added complexity of having 3 branches in every function. Why not just use the generic case everywhere?

The main reason is to avoid numerical noise (I want the cosine of a real number to be a real number). You have

sage: a = CIF(2,0)
sage: I = CIF.gen()
sage: ((I*a).exp() + (-I*a).exp()) / 2
-0.4161468365471424? + 0.?e-16*I

It should be fine to remove the imaginary branch, but then it would be less precise

sage: a = CIF(0,2)
sage: b1 = CIF(0, a.imag().cosh())
sage: b1.diameter()
1.33393183261575e-16
sage: b2 = ((I*a).exp() + (-I*a).exp()) / 2
sage: b2.diameter()
2.36079803558624e-16

Vincent

comment:11 Changed 5 years ago by jdemeyer

  • Milestone changed from sage-6.4 to sage-6.6
  • Status changed from needs_review to needs_work

For the sine and cosine, you can use the formula at the end of http://en.wikipedia.org/wiki/Trigonometric_functions#Relationship_to_exponential_function_and_complex_numbers which computes the real and imaginary part separately. Assuming that sinh(0) = 0 exactly, this formula will not introduce a fake imaginary part.

comment:12 Changed 5 years ago by git

  • Commit changed from dd4cf4e61cbc44485985e1f74de201c87ad06bc9 to 3f8cd4178d279dbb6e5ad87c3f4a1c9f5ac75548

Branch pushed to git repo; I updated commit sha1. This was a forced push. New commits:

3f8cd41Trac 17285: CIF trigonometric functions

comment:13 Changed 5 years ago by vdelecroix

  • Status changed from needs_work to needs_review

New implementation using other trigonometric identities (and much faster).

For tan and tanh this is certainly not the best possible.

Vincent

comment:14 Changed 5 years ago by jdemeyer

  • Reviewers set to Jeroen Demeyer
  • Status changed from needs_review to needs_work

You're not actually supposed to do this:

from sage.ext.interrupt cimport sig_on, sig_off

The official way is still include "sage/ext/interrupt.pxi" (and this will even become mandatory after #18027)

comment:15 Changed 5 years ago by jdemeyer

I like the new code. However, a few comments:

  • move the # cos(x + yi) = ... comments inside the docstring as an ALGORITHM: block.
  • move the mpfi_init2() and mpfi_clear() outside the sig_on()/sig_off() (it's not really wrong to put them inside, but it's certainly safer outside)

comment:16 Changed 5 years ago by git

  • Commit changed from 3f8cd4178d279dbb6e5ad87c3f4a1c9f5ac75548 to afccfe62eac92dd2cea900f3bac44600e633aa7d

Branch pushed to git repo; I updated commit sha1. New commits:

afccfe6Trac 17285: signal handling + doc

comment:17 Changed 5 years ago by vdelecroix

  • Status changed from needs_work to needs_review

comment:18 Changed 5 years ago by jdemeyer

  • Status changed from needs_review to positive_review

comment:19 follow-up: Changed 5 years ago by vdelecroix

Thanks for the review.

Just one related question. With the current code, if a KeyboardInterrupt happens in between the sig_on/sig_off block, then the code mpfi_clear(tmp) will never get called and hence the memory allocated to tmp will not be freed. Is there a clean solution for that?

comment:20 in reply to: ↑ 19 Changed 5 years ago by jdemeyer

Replying to vdelecroix:

Is there a clean solution for that?

try:/finally:.

comment:21 Changed 5 years ago by vdelecroix

  • Description modified (diff)

comment:22 Changed 5 years ago by vbraun

  • Branch changed from u/vdelecroix/17285 to afccfe62eac92dd2cea900f3bac44600e633aa7d
  • Resolution set to fixed
  • Status changed from positive_review to closed
Note: See TracTickets for help on using tickets.