#20438: Wrapper_cdf.call_c() drops imaginary part
-------------------------------------------------+-------------------------
Reporter: cswiercz | Owner:
Type: defect | Status: new
Priority: major | Milestone: sage-7.2
Component: cython | Resolution:
Keywords: fast_callable, interpreters, | Merged in:
wrapper_cdf | Reviewers:
Authors: cswiercz | Work issues:
Report Upstream: N/A | Commit:
Branch: | Stopgaps:
Dependencies: |
-------------------------------------------------+-------------------------
Comment (by cswiercz):
Of course! Below is the Cython input.
{{{
# input
# is there no Sage %%cython cell magic?
cython("""
from sage.ext.interpreters.wrapper_cdf cimport Wrapper_cdf
# cannot do below: it redefines double_complex from the
# auto-generated headers created by ext.interpreters.py. generates
# cython error:
#
# "Cannot assign type 'double complex [2]' to 'double_complex *'"
#
#ctypedef complex double_complex "double complex"
# Instead, use the (incorrect) definition as given in the auto-gen
# wrappers. This ctypedefs "double_complex" as a "double"! That is,
# all complex parts of calculations are lost when invoking
# Wrapper_cdf.call_c()!
#
# The following is the self-proclaimed "hack" in the auto-gened file
# ext/interpreters/wrapper_cdf.pyx
cdef extern from "complex.h":
ctypedef double double_complex "double complex" # cannot just change
"double" to "complex"
cdef double creal(double_complex)
cdef double cimag(double_complex)
cdef double_complex _Complex_I
cdef inline complex call_c(Wrapper_cdf f, complex x, complex y):
# Quickly evaluates a fast_callable function `f` on domain `CDF`
# (a Wrapper_cdf tupe) with complex arguments `x`, and `y`.
cdef double_complex[2] input
cdef double_complex* output = [0]
cdef complex value
input[0] = <double_complex>x
input[1] = <double_complex>y
f.call_c(input, output)
value = <complex>(output[0])
return value
from sage.all import QQ, CDF, fast_callable
R = QQ['x,y']
x,y = R.gens()
f = y**3 - 2*x**3*y + x**7
f_fast = fast_callable(f, vars=[x,y], domain=CDF)
cdef complex a = 1. + 1.j
cdef complex b = 1. - 2.j
cdef complex value_sage = f_fast(a,b)
cdef complex value_c = call_c(f_fast,a,b)
print 'sage value: ', value_sage
print 'call_c value:', value_c
""")
}}}
And the corresponding output of the last two print statements:
{{{
sage value: (-7-18j)
call_c value: (-7+0j)
}}}
--
Ticket URL: <http://trac.sagemath.org/ticket/20438#comment:2>
Sage <http://www.sagemath.org>
Sage: Creating a Viable Open Source Alternative to Magma, Maple, Mathematica,
and MATLAB
--
You received this message because you are subscribed to the Google Groups
"sage-trac" group.
To unsubscribe from this group and stop receiving emails from it, send an email
to [email protected].
To post to this group, send email to [email protected].
Visit this group at https://groups.google.com/group/sage-trac.
For more options, visit https://groups.google.com/d/optout.
