expose sincos?

[mlubin]

GNU libm and openlibm both implement sincos:

void sincos(double x, double *sin, double *cos);


DESCRIPTION
       Several  applications need sine and cosine of the same angle x.  This function
       computes both at the same time, and stores the results in *sin and *cos.

LLVM also has a sincos intrinsic.

Could we expose this? I observed a 30% speedup by calling sincos over sin and cos separately when using the openlibm implementation. Ideally this would be a compiler optimization, but since we're not using the LLVM intrinsics for sin and cos, I'm guess that we won't easily get LLVM to do this for us for free.
In DualNumbers, this would be very useful since computing either sin or cos of a Dual requires both values, which in turn would speed up derivative evaluation in JuMP.

[mlubin]

I would submit a PR for this, but I don't know how to wrap the libm signature of sincos without allocating a temporary array for the output.

[simonbyrne]

Some previous discussion on julia-users:
https://groups.google.com/d/topic/julia-users/368fpIGYRrY/discussion

I think the easiest solution might be simply to implement these in Julia: we already have the sin/cos "kernel" functions, as they are used in sinpi/cospi, so all that's really needed is the range-reduction code (which would be useful for #8799 anyway). I have made a start here:
https://github.com/simonbyrne/libm.jl/blob/master/src/rempio2.jl
It still requires an implementation of Payne-Hanek reduction (required for arguments > 1.6e6) but until that gets done, we could just call the existing C code.

[mlubin]

@simonbyrne, could you give a hint on how to put all the pieces together to test this out if I don't care about arguments > 1.6e6?

[simonster]

We should also use this for cis, which it looks like currently calls sin and cos separately.

[simonbyrne]

@mlubin I had something in mind like this:
master...simonbyrne:sincos

It still does seem to be a bit slower, but I haven't really done any profiling to find out where the problem is.

UPDATE: the overhead is the array allocation for ccall, even when that branch is never called. If it is commented out, I see a small (10% or so) speed up over separate sin/cos calls.

[mlubin]

Cool! Is there any extra speedup from combining the

    z = x.hi*x.hi
    w = z*z

computations across cos_kernel and sin_kernel?

[simonbyrne]

From what I can tell, if you @inline the kernel functions, then LLVM is able to combine them itself. This does give a nice further speedup.

[vtjnash]

I would submit a PR for this, but I don't know how to wrap the libm signature of sincos without allocating a temporary array for the output.

with the ccall work that's being merged soon, you can replace 1-value Arrays with Ref objects. eventually, those should even get allocated on the stack too.

[zsoerenm]

This can also give a tremendous speedup for exp(1im * arg)
as it is the same as complex(cos(arg), sin(arg))
Further details are discussed here: https://discourse.julialang.org/t/why-is-this-julia-code-considerably-slower-than-matlab/2365/45

[yuyichao]

This requires better optimization before it can be implemented.

[StefanKarpinski]

I'm unclear on why we can't just expose sincos and do whatever other stuff is suggested in this thread.

[yuyichao]

What's the "other stuff is suggested"?

[StefanKarpinski]

I.e. why can't we expose a sincos function before optimizing it fully? Then we can work on optimization at any point we feel like it?

[yuyichao]

We can define sincos, which we kind of do already and it's called cis. We can't expose sincos from libm since that'll be slower.

[StefanKarpinski]

Well, cis kind of computes sincos if you understand complex exponentials. I suppose if we allowed destructuring of complex numbers then one could just write this: c, s = cis(x). Actually, that would be fairly handy for rationals as well.

[simonbyrne]

meh: destructuring complex numbers seems a little too magical. I don't think it's too hard to do reim(cis(x)). The other issue is that it gives them in the opposite order.

I would be in favour of a sincos function (as well as sincospi and possibly a sincosd as well).

[Datseris]

Hi,

sorry for posting on such an old issue, but I also need the function sincos for a library I am using. However, even though this issue thread seems to report that the function is accessible (and in fact it is in julia/exports.jl), for some reason I don't have access to it:

   _       _ _(_)_     |  A fresh approach to technical computing
  (_)     | (_) (_)    |  Documentation: https://docs.julialang.org
   _ _   _| |_  __ _   |  Type "?help" for help.
  | | | | | | |/ _` |  |
  | | |_| | | | (_| |  |  Version 0.6.2 (2017-12-13 18:08 UTC)
 _/ |\__'_|_|_|\__'_|  |  Official http://julialang.org/ release
|__/                   |  x86_64-w64-mingw32

julia> sincos
ERROR: UndefVarError: sincos not defined

help?> sincos
search:

Couldn't find sincos
Perhaps you meant sinc, sin, sind, sinh, sinpi, sizeof, cos, acos, minabs, info or asin
  No documentation found.

  Binding sincos does not exist.

I am on Windows 10 64-bit. Can you tell me where I can find the function?

[fredrikekre]

Can you tell me where I can find the function?

#21589 was merged post 0.6. It's on julia master.