Currently it yields the integer 0.
Should there be an sgn() of an imag(w) or a complex(w).
Should sgn() of a floating point type be integer or a real(w)?
Currently
sgn(+0.0) == sgn(-0.0)
which is arguably wrong.
I do not think it is all that urgent but it might mean that sgn() is another @unstable candidate.
RE: sgn(0.0) == sgn(-0.0)
I would tend to disagree this is wrong, both 0.0 and -0.0 are zero and the sign function is defined to be zero at zero, hence I think it's correct sgn(0.0) == sgn(-0.0) == 0.0
I guess we could have
sgn(0.0) -> 0.0
and
sgn(-0.0) -> -0.0
but feels a bit overengineering maybe? Also since 0.0 == -0.0 returns true, this would not change the original equality (which I think should hold)
continue: sorry I prematurely pressed enter in the previous message
Should there be an sgn() of an imag (w) or a complex (w).
I am not sure if it's useful or not, but in principle the sign function is generally extended to complex numbers using the definition sgn(z) = 0 when z=0 and sgn(z) = z/|z| otherwise.
should sgn () of a floating point type be integer or a real (w)? what should sgn(NaN) be?
I believe sgn of real(w) should be real(w) and my motivation for this actually stems from your second question. I believe we want sgn(NaN) to return NaN. For the simple fact that we want NaN to propagate, e.g. I could have sgn(some_expression(x)) where some_expression(x) produces NaN because of some issues in my code.
If we want sgn(NaN) to return NaN then this implies that sgn should return real (otherwise if it returned an integer sgn(NaN) should error or return some arbitrary integer, both of which sound a bit meh).
Thanks for your comment Luca. It has made me realise my original post was a bit sloppy in its wording..
Your commented that
sgn(+0.0) -> +0.0
sgn(-0.0) -> -0.0
feels a bit over engineered maybe.
The IEEE 754 (ISO 60559) standard is built on top of the fact that +0.0 and -0.0 are ... (at the risk of being way too brief) different in assignment while behaving as if they were the same in a comparison. So, I think Chapel has to live with what can be called over engineering.
That said, your reading of my words showed that my original text was imprecise at best, or wrong at worst. I should have said that
sgn(+0.0) should return +0.0
sgn(-0.0) should return -0.0
So, while the two results are not the same bit pattern in memory because they differ in their negative field, they compare as equal computationally, i.e.
sgn(+0.0) == sgn(-0.0)
because
+0.0 == -0.0
Both values co-erce to the same integral value
.
The 1-bit wide negative field is called the sign bit in some fields.
Thanks for your second post. I wanted to canvas others' ideas and explanations and expectations, ... especially as I am suggesting a change of return type on Chapel's existing
inline proc sgn(x : real(?w))
and adding overloads for imag(w) and complex(w) for completeness.
Like you, I am a big fan of NaN propogation.
We might see what others say eh before I waffle on any more.
We already have proc signbit(x : real(?w)) that is unstable. In principle, we could make that one return real as described here (and in issue signbit() - simplification · Issue #24578 · chapel-lang/chapel · GitHub) to get the sign of a zero and to propagate NaNs and we can do that after the 2.0 release (which is coming up fast).
But, this thread is proposing that proc sgn(x : real(?w)) be unstable and/or return a real.
So my question is: would it be sufficient to say something along the lines of "use proc signbit if you want NaN propagation" ? If it is sufficient, we don't need to change proc sgn in this way. If it isn't sufficient, why not?
I have never seen a case where I needed to worry about NaN propagation with signbit(), mainly because it is either irrelevant totally or I have covered it myself with separate logic.
This thread is proposing marking sgn(x : real(?w)) as unstable with a view to changing its type to return a real(w).
How do I get that cool blue text like you have just done?
By changing its type to real(w), you also get to return a zero result with the correct (and IEEE 754 compliant) signed zero (which I often need and which I would hesitate to suggest others might also need - even if they do not know it).
As per the following
inline proc sgn(x : imag(?w)) : real(w) ...
inline proc sgn(x : complex(?w)) : real(w/2)...
the other two extra proposed routines also return respectively real(w) and real(w/2) numbers.
We already have
proc signbit(x : real(?w))that is unstable. In principle, we could make that one returnrealas described here (and in issue signbit() - simplification · Issue #24578 · chapel-lang/chapel · GitHub) to get the sign of a zero and to propagate NaNs and we can do that after the 2.0 release (which is coming up fast).
I would find this a bit confusing, since it seems like it mixes signbit and sgn. First, I would expect the function signbit to return a bool based on its name. Having signbit to return real (which propagates NaN and handles 0's signs etc.) while sgn returns an integer sounds a bit strange. Also, would signbit in that case return 1.0 and 0.0?
To me, it would make more sense to have signbit return only the signbit of the number as bool and sgn behave more like the mathematical sign function, behaving as damian described (return an output of the same type, propagate nan etc.). An added advantage of this is that I could do things like if signbit(x) then
The following implementation could work (quick prototype, there's probably room for more low-level optimization. The procs should probably be inlined and it might be worth rewriting with less if statements (not sure would need to benchmark))
proc my_sgn(x: ?T): T where isComplex(T) {
if x == 0 then return x;
else return x / abs(x);
}
proc my_sgn(x: imag(?w)): imag(w) {
return my_sgn(x: real(w)): imag(w);
}
proc my_sgn(x: ?T): T where isIntegral(T) {
if x == (0: T) then return x;
else if x < 0 then return (-1: T);
else return (1: T);
}
proc my_sgn(x: ?T): T where isReal(T) {
if x == (0: T) || isNan(x) then return x;
else if x < 0 then return (-1: T);
else return (1: T);
}
in principle, the first proc for isComplex would work for all the cases (using the isNumeric clause), but that would probably cost some performance, so would make sense to have more functions for real types.
I think they are markdown links typed as [text to display](link), for example
[`proc signbit(x : real(?w))`](https://chapel-lang.org/docs/main/modules/standard/Math.html#Math.signbit)
would be displayed as proc signbit(x : real(?w))
@damianmoz and @lucaferranti - thanks very much for your feedback here. I've gone ahead and created PR #24582 to make these unstable. Based on some offline discussion, I'm expecting to get this merged just in time for the Chapel 2.0 release. Because they will be unstable, we'll be able to make adjustments to them in the future.
For the longer-term question of how to adjust them, I've created issue #24583 to discuss further.
Thanks Michael.
Thank you too Luca. Anytime you want a discussion about signed zero, a bedrock feature of the IEEE 754 standard, send me an email. Sadly, we have not managed to communicate the importance of signed zero to end users, I was thrilled that you were worried about NaN propagation.
Unlike #19165 on Github, this posting focused on IEEE 754 behavior with signed zeros and NaNs, two issues which came up recently on two forums with implications and ramifications broader than just the Chapel community. I thought it would be useful and opportune to leverage what was happening in those two forums because the amount of work to be done in Chapel seemed small and I could do most of it so it was not like I was overloading anybody else.
By the way, an IEEE 754 binary floating point number is defined as
x = (-1)^n * 2^e * 1.f
where n is the negative field, e is the exponent field and 1.f is the significand although only the bits if f are stored in memory. Some people call the tiny integer n the sign bit. So, signbit() is the way to capture an integral value, not a bool. More correctly, an unsigned integral value. Some implementations follow which have a lot less overhead that your earlier examples because they exploit Brad's transmute() primitive - thanks again Brad!!.
inline proc sgn(x : real(?w))
{
type Real = real(w);
type Uint = uint(w);
param _one = (1:Real).transmute(Uint); // get the bit pattern for 1.0:real(w)
param _neg = 1:Uint << (w - 1); // form a bit mask to extract the negative field
const _x = x.transmute(Uint); // get the bit pattern of the argument
return if x != x || (_x << 1) == 0 then // when x is a NaN or x == +/- 0, just return the original
x
else // grab the signbit from the top bit of x and stuff it into the top bit of +1.0
(_one | (_x & _neg)).transmute(Real);
}
inline proc sgn(y : complex(?w))
{
param _w = w / 2;
type Real = real(_w);
type Uint = uint(_w);
param _one = (1:Real).transmute(Uint);
param _neg = 1:Uint << (_w - 1);
const x = y.re;
const _x = x.transmute(Uint);
return if x != x then // i.e. when x is a NaN, just return the original
x
else if (_x << 1) == 0 then // when x == +/- 0, look at the imaginary component
sgn(y.im)
else // stuff the signbit() of x into the top bit of +1 as seen above
(_one | (_x & _neg)).transmute(Real);
}
inline proc sgn(y : imag(?w))
{
return sgn(y:real(w));
}
Sorry, the param versions are missing but hopefully you get the drift.
The C version of the first of these compiles into 14 X86-64 instructions, none of which involve arithmetic. So, pretty lean and mean. I have not checked the Chapel.
Hehe, I am a floating point geek, I bought myself the Handbook of floating-point arithmetic as christmas gift a couple of years ago 
always happy to chat about the topic!
Your implementation of sgn for real looks great, I think it's the way to go! I was thinking something on the same lines with bit operations and was looking for something like transmute but could not find it, thanks for pointing that out!
I am a bit less sure about the proposed sgn(z: complex(?w)): real(w). What's the intended semantics of the function? It seems to look at sgn for the real part first and then at the imaginary part if the real part is zero. What is the rationale behind this?
Alternatively, what would you think about defining sgn(z: complex(?w)): complex(w) to basically compute z/abs(z) (returning z when z=0, not sure what would make sense to return if the real and/or imaginary part is +/- inf). This is generally how the sign function is defined for complex numbers (similar for imag(?w)). This extends the invariant x = |x| * sgn(x) that holds for real numbers to complex numbers as well (well modulo rounding error of course).
RE: signbit
I agree that from the equation (-1)^n * 2^e * 1.f then n should be an integer (constrainted to be either 0 or 1). Having n as bool then to be rigorous the equation should probably be adjusted to be (-1) ^ [n] * 2 ^ e * 1.f, where [...] is iverson bracket. To me this is a distinction between a data-structure and its semantics. In this sense I would say a floating-point is defined as a triplet (n, e, f) and represents a real number as you wrote (plus all peculiarities, like subnormal numbers, 0, exponent bias, nan, inf, etc.). This is also a distinction made in the standard (2019 revision, table 3.1, although the layer 3 triplet in the standard is a bit higher level than mine here, I guess my triplet would be a layer 3.5) (I know I am getting a bit pedantic here, sorry, told you I like chatting about floating point arithmetic )
The standard (section 5.7.2, 2019 revision) says that an implementation shall provide a function isSignMinus that returns a bool and is true if and only if the input has negative sign (applying also to NaN). So to me the signbit function is this (this is also the choice in C++ and julia). Now I don't have strong opinions of what this function should be called (e.g. rust calls it is_sign_negative), but if sgn(x: real) behaves as you outlined (which I think it should!) then I am not sure we need another function returning an integral for the sign bit (but maybe I'm missing some valid use-cases), but we do need one returning a bool, and the name signbit would make sense to me.
Before I reply to your most recent post, once Chapel supports something like C's static const, you can get even tighter code.
proc sign(x : real(?w))
{
type Real = real(w);
type Uint = uint(w);
param one = 1:Real;
param signedOnes : [0..1:Uint] Real = [ one, -one, ];
const _x = x.transmute(Uint);
return if x != x || (_x << 1) == 0 then // x is a NaN or x == +/- 0
x
else
signedOnes[_x >> (w - 1)];
}
Hi Luca,
Brad implemented a transmute method last year. I had done a crude version (which only worked with const, not param data) a few years prior. I will send you some rough doco over the weekend (unless Brad has something already). Given a real(w) variable v, then
const v_as_bits = v.transmute(uint(w));
will return an unsigned integral type which is the same bit pattern as that real(w) number v. Also
const v_as_real = v_as_bits.transmute(real(w));
assert(v == v_as_real);
assert(v == (v.transmute(uint(w)).transmute(real(w))));
Your definition of sgn() for a complex is the original signum() definition for this functionality which I think has been around since at least the early 1800s if not a little bit earlier. The name sgn() did not originate for this functionality until the very early 1900s. Kronecker in about 1887 used [k] to represent signum(k) for integers and his definition used the form with which we are familiar today. I hate the name sgn() but that is another story. I think we can blame an Italian, Peano for that name but I am not so sure.
The implementation I provided for sgn(z : complex(w)) is based on a formula I have used for years in practice in my programming of engineering, physics and
geophysics applications. I believe Maple uses this same definition and it is fast for a complex number. Matlab however uses the classical definition that you quote and is slow. What is your application? I find the classical definition impractical in my work. And besides, I can write
inline proc signum(z) // totally generic
{
return if z == 0 then z else z / abs(z);
}
if I really want the classical routine. A bit slow though.
The negative field n is an integer. It is never a bool. The equation as I wrote it is correct. That n is only ever 0 or 1.
And yes, the standard says there shall be a routine isSignMinus. Chapel is your friend. If you really want this routine, and I have never needed a bool, provide your program with
inline proc isSignMinus(x : real(?w)) do return signbit(x) != 0;
The routine signbit() is a better and faster building block while isSignMinus() is not well suited to being the building block.
The sgn() function never returns the value of the negative field, i.e. n, or what most people call the sign bit. As I often want the value n I need a routine like signbit() which returns that integer. If you want a bool, use isSignMinus as described above (which is also what Swift uses).
Hi Damian —
Regarding:
There is an early draft of this feature going in as part of Chapel 2.0.
However, it doesn't quite get us where we'd want for this use-case because in its current form, it only supports a singleton value on the first locale to encounter the definition. A next step in the effort is to extend it to support per-locale static values at which point it should be ready for math routines like this one (where I'm imagining you'd want to apply it to signedOnes in the absence of actualy param arrays, is that right?).
-Brad
Thanks for the update Brad.
Yes, whatever the feature is which provides the equivalent of something that is not initialized every time the routine is called, would be used for signesOnes.
The most common case would be an array of probably polynomial coefficients, which could have real(w), imag(w) or even complex(w) type. Less common might be an array of tuples or records. All the data would be known at compile time.
Is this new feature evaluated at compile time or first time the routine is used. Or is it smart enough to realise what it can do at compile time to avoid any evaluation? Ta.
@damianmoz : The current draft of static local variables initializes the variable the first time its declaration is encountered. The assumption is that anything that can be evaluated at compile-time would be expressed as a param (and nothing about this effort affects the current bounds between what the param machinery can and cannot handle).
-Brad
Got it. I might just stick with a private const defined outside the routines that need these little arrays for now as it seems like this new feature has a run-time overhead and might affect the optimization of the routine.
There is no rush with a param of a structured type or an array. But I think Chapel needs it eventually. Everything from a fundamental mathematical library to major HPC programs use lookup tables which need to have effectively a zero run-time overhead. Sometimes you can fudge it with params like c0, c1, c2, ..., cN`. But in other places you need to be able to index an array.
Hi Damian —
I think that that's a good alternative in the short-term. Part of the motivation for the feature is that that workaround doesn't work in the context of a generic routine where you'd need a copy of the local variable per instantiation of the routine. But for non-generic routines in multi-locale settings, the module-scope private const approach is currently preferable due to the per-locale replication it gets you. I can't speak to the relative overheads of the two approaches or their impact on optimizations at present, but that'll definitely be something to keep an eye on.
-Brad