What we need

[simonbyrne]

The following functions are defined in the C 2011 standard, under math.h. We probably don't need all of these as some (such as trunc or round) are better provided by LLVM intrinsics.

Trig functions

• acos
• asin
• atan
• atan2
• sin
• cos
• tan

Hyperbolic functions

• acosh
• asinh
• atanh
• cosh
• sinh
• tanh

Exponential and logarithmic

• exp
• exp2
• exp10
• expm1
• ilogb
• ldexp
• scalbn
• log
• log10
• log1p
• log2
• logb
• frexp
• modf

Power and absolute value

• cbrt
• pow
• hypot
• erf
• erfc
• lgamma
• tgamma
• fabs: use intrinsic
• sqrt: use intrinsic

Nearest integer

• ceil
• floor
• round
• trunc
• rint
• nearbyint

Remainder functions

• fmod
• remainder
• remquo

Manipulation

• copysign: use intrinsic
• nan
• nextafter
• nexttoward

Maximum, minimum

• fmin
• fdim
• fmax
• fma: use intrinsic
  • software fma (corner case where system libm fma is not sufficient)

[ararslan]

Thanks for putting this together! This will certainly be a fun one. 😄

Any idea how much of your Libm.jl code you're planning to bring over here, if any?

[quinnj]

A few of these are already implemented in Base in native Julia, right? We should probably port those over here or at least check them off the list here.

[yuyichao]

Do we not need bessel functions anymore?

[musm]

How do I access the intrinsics?

[yuyichao]

Base.llvmcall

[simonbyrne]

I should say, this wasn't intended to be a definitive list:

• we don't have to keep the names, i.e. it probably makes sense to keep fabs as abs, and tgamma as gamma
• there's probably no point in including the ones that use intrinsics in here (that just seems like wasted effort)
• we may want to move some of these functions to SpecialFunctions.jl
• conversely, there may be other functions we would like to include (e.g. those in IEEE-754 recommended functions JuliaLang/julia#6148).

[simonbyrne]

Do we not need bessel functions anymore?

Personally, I think these might be a good candidate for SpecialFunctions.jl

[simonbyrne]

Any idea how much of your Libm.jl code you're planning to bring over here, if any?

I'll have to look: a lot of it is quite old, but if you see anything you would like, feel free to copy it over.

[musm]

Base is using a julia implementation of hypot https://github.com/JuliaLang/julia/blob/2d24eda52727c6cbccec935d689d84d710366f5a/base/math.jl#L194

[ararslan]

Personally, I think these might be a good candidate for SpecialFunctions.jl

They're already there. 😄

[ararslan]

I think the error and gamma functions would make more sense in SpecialFunctions.jl.

[oxinabox]

Base is using a julia implementation of max and min
https://github.com/JuliaLang/julia/blob/2d24eda52727c6cbccec935d689d84d710366f5a/base/math.jl#L234

so fmax and fmin can probably be struck.

[oxinabox]

Julia is using the LLVM intrinsic remainder (frem) already.

It is Core.Intrinsics.rem_float
referenced at https://github.com/JuliaLang/julia/blob/2d24eda52727c6cbccec935d689d84d710366f5a/base/float.jl#L254

Julia also has its nan's defined. So the C nan function/constant is not needed
https://github.com/JuliaLang/julia/blob/2d24eda52727c6cbccec935d689d84d710366f5a/base/float.jl#L30

[simonbyrne]

Actually fmin and fmax are different from Julia's min/max due to how they handle NaNs, see JuliaLang/julia#7866

[musm]

on the checkbox fma: there is no llvm intrinsic if the cpu doesn't support fma, we will actually need a a fma function too

[vchuravy]

The intrinsic still exist, but it will not be lowered to a single hardware instruction. http://llvm.org/docs/LangRef.html#llvm-fma-intrinsic
and we have it available at Core.Intrinsics.fma_float

[musm]

My bad.
Currently we either use hardware fma if available or the openlibm fma function if the hardware fma is not available, and not the llvm intrinsic that you linked, if I am not mistaken.

[vchuravy]

fma_llvm uses the Core.Intrinsics version to lower to the hardware fma see https://github.com/JuliaLang/julia/blob/aebf3ac00eff4841a034a2f8c4e2d0263da3669d/base/floatfuncs.jl#L255-L269

[musm]

Right, in that file it selects the software fma provided by openlibm if the intrinsics are not accurate enough. At the end of the day we need a software fma in pure julia to replace that, which is what I meant to say.

[simonbyrne]

Actually, it's a little more complicated than that. The LLVM fma intrinsic will fall back to the LLVM-linked libm (which is usually the system libm). The problem is that on old glibc versions, this clobbered the rounding mode which caused all sorts of havoc (see JuliaLang/julia#9847). So instead we now call the openlibm fma function if we don't think there's a native fma (which is what that line you link to is trying to test).

[musm]

Aren't
log2 and logb the same function
similarly
scalbn and ldexp
unless we want to support system's where the float base is not 2?

[musm]

Is there an llvm intrinsic that does round nearest away (i.e. 2.5 ->3.0 and -2.5 -> -3.0)? rint_llvm rounds nearest towards.

Currently, I use

rint{T<:FloatTypes}(x::T) = ifelse(x < 0, x - T(0.5), x + T(0.5))

But it's not as good as a single vroundsd instruction using round

[musm]

Looking at: http://support.amd.com/TechDocs/26568.pdf there doesn't seem to be an instruction that does round away from zero

Edit: upon closer inspection, there isn't actually a need in any of the algorithms to have the rounding behave this way. As a result, i'm transitioning to use round in the code.

[oxinabox]

Should tgamma, and lgamma be struck out,
as they are currently in Base, in pure julia?
(and will, (probably) at somepoint, move to SpecialFuns.jl)

And erf and erfc are currently in openspecfun,
which suggests it is in scope of SpecialFuns.jl
to re-implement them in pure julia,
rather than here.
So maybe those can be knocked out too?

[StefanKarpinski]

Dang. You guys have been really making moves here.

[musm]

@StefanKarpinski this is largely irrelevant now for Amal, and only applies to https://github.com/JuliaMath/Sleef.jl

I should close this and open another one relevant to Amal.

Sleef has important limitations in regards to performance and that was the reason I decided to start from scratch.