41/42/...
Challenge: Create a program that will return the number of significant digits.
When run on a 41, it must return 10. The *same* program run on a 42S must return 12, and on Free42 decimal, it must return 25.
Cheers, Werner
01>LBL "DIGITS"
02 1
03 1.999
04 STO 00
05>LBL 05
06 Rv
07 10
08 *
09 ENTER
10 ENTER
11 ENTER
12 1
13 +
14 
15 X=0?
16 GTO 19
17 ISG 00
18 GTO 05
19>LBL 19
20 RCL 00
21 IP
22 RTN
How about this?
Lyuka
My favorite one:
LBL "T"
1
3
1/x
3
*

LOG
CHS
These have been tested on the HP42S, Free42 Decimal and on the HP15C:
00 { 17Byte Prgm }
01>LBL "SD"
02 6
03 ENTER
04 1/X
05 1/X
06 
07 1/X
08 LOG
09 1
10 +
11 .END.
001 f LBL A
002 6
003 ENTER
004 1/x
005 1/x
006 
007 1/x
008 g LOG
009 1
010 
011 g RTN
Cheers,
Gerson.
Oops! The step 010 in the HP15C program is
010 +
That's been a typo, obviously :)
Hi, Werner:
This doesn't meet the requirements of your challenge as it's not meant for the 41/42/... as asked but for the HP71B instead, however I thought my approach might result novel.
For the HP71B you don't actually need a program, just execute this from the command line:
LGT(1NEIGHBOR(1,0))
12
This mathematical expression will work in any other calculator or computer provided you express it in its native paradigm (RPN, RPL, ...) and have available (or can inmplement) the two numeric functions required, namely LGT (base10 logarithm) and NEIGHBOR (machine successor).
Best regards from V.
Same logic On a 50G :
« 1 3 INV 3 *  LOG NEG »
or
« 3 INV >STR TAIL SIZE »
Edited: 17 Jan 2012, 12:05 p.m.
Nice.
Insert >NUM to make this work in exact mode, too:
« 1 3 INV 3 *  LOG NEG >NUM »
« 3 INV >NUM >STR TAIL SIZE »
If you EVAL the first, instead of doing >NUM, it returns +Infinity, as it should.
On ND1, it returns Infinity, too, even though there's no exact mode.
Intriguingly, 0.3333... * 3 is 1 with IEEE754!
Quote:
Challenge: Create a program that will return the number of significant digits.
The significant digits
of the machine, that is. When I first read this I assumed that the number of significant digits
of X should be returned.
In other words: Enter 0.001234 or 1234 or 1234000 and get 4, while Pi would return 10 or 12, depending on the machine. Finally, 0 should return 0 ...or maybe 1?
Any solutions for this one ?)
Dieter
For the HP 50g:
<< MANT >STR SIZE 1  >>
Gerson.
On the WP 34S you need to add an ENTER^ after the first number or replace it by # 1 (from the CONST catalogue). It returns 16 in single precision and 34 in double precision which hits the nail on the head.
Quote:
On ND1, it returns Infinity, too, even though there's no exact mode.
Intriguingly, 0.3333... * 3 is 1 with IEEE754!
Are you doing binary arithmetic in ND1? With decimal arithmetic this shouldn't happen.
In single precision this can be done on the 34S, too:
ALL 00
CL[alpha]
MANT
[alpha]RC# X
[alpha]LENG
DEC X
For double precision the alpha register is too short.
Edit: Thinking again this does not work :( [alpha]RC# X does format to the display size of 12 digits, not the full precision. Back to the drawing board!
Edit again: A loop should do:
MANT
0
x[<>] Y
LBL 00
x[!=]0?
SKIP 02
DROP
RTN
INC Y
FP
SDL 01
GTO 00
Edited: 18 Jan 2012, 6:29 a.m. after one or more responses were posted
This works:
LBL A
# 1
9
1/x
MANT
INC Y
SDL 01
FP?
BACK 003
DROP
END
 Pauli
Marcus,
You mean BCD, right? ("Binary arithmetic" sounds like arithmetic with binary numbers...) Sadly, ND1 doesn't have a BCD option yet.
Reals are IEEE754 doubles in ND1. And the RPL programs above imply reals. (Unlike on the 50g in exact mode, where they imply BigInts. ND1 does have BigInt support, but I chose to not adopt the HP 50g convention (in exact mode only) to regard a decimal typed without a "dot" as a BigInt. Instead, if you want a BigInt, you have to use R>I (or toBig, in "modern" syntax).)
In the latest release BigFloats are supported, too.
So,
« 3 toBigF INV >NUM >STR TAIL SIZE »
will actually return whatever the BigFloat precision is set to (40 digits, by default).
To confirm that "1/3*3" is "1" in IEEE754, just type this in gdb:
p 1.0/3.0*3.0
Works for any integervalued double I tried.
Not bad for a format, for which 1.12 is 1.1200000000000001.
(Yes... I'm still yearning for including BCD at some point.)
EDIT: I educated myself and do understand now that with "binary arithmetic" you meant (and correctly termed) what to me is floatingpoint arithmetic. So, yes, ND1 is using binary arithmetic exclusively, for its Real data type.
Edited: 18 Jan 2012, 3:33 p.m.
Quote:
Any solutions for this one ?)
Here's my attempt at the HP42S version. It will work for 0.001234 or 1234 or 1234000 but not for pi (because it's limited to 11 significant digits only)*. 0 will return 1 (because I think it should— but I may be wrong). On the HP41 (not tested), all 1E10 constants should be replaced with 1E08.
00 { 33Byte Prgm }
01>LBL "SD"
02 XEQ "MANT"
03 1
04 X<>Y
05>LBL 00
06 FP
07 X=0?
08 GTO 01
09 10
10 ×
11 1
12 STO+ ST Z
13 Rv
14 GTO 00
15>LBL 01
16 Rv
17 .END.
00 { 37Byte Prgm }
01>LBL "MANT"
02 X=0?
03 RTN
04 ABS
05 LOG
06 FP
07 10^X
08 1E10
09 ×
10 0,5
11 +
12 IP
13 1E10
14 ÷
15 1
16 X>Y?
17 10^X
18 ×
19 END
Gerson.

* 12345678910 XEQ SD currently returns 11 instead of 10 because MANT returns 1.2345678909. Too bad the HP42S doesn't have a MANT command like the WP 34S—
Edited: 17 Jan 2012, 10:46 p.m.
It needs one small further edit :)
 Pauli
Almost what I had in mind..
01*LBL"SD"
02 .6
03 ENTER
04 1/X
05 1/X
06 
07 LOG
08 CHS
09 END
Werner
use this MANT:
*LBL"MANT"
X=0?
RTN
ENTER
ABS
LOG
INT
LASTX
X<Y?
DSE Y
LBL 00
RDN
10^X
/
RTN
It only ever scales by a power of 10, and will thus never lose a digit.
Cheers, Werner
To get this working on the 30b you need to input INPUT twice (or should I have said 'enter ENTER twice'? ;)). The result is 12.
Quote:
Here's my attempt at the HP42S version.
This one uses another approach, here for the 41/42.
01 LBL"SDIG"
02 ,1 ; dummy ISG control number
03 X<>Y
04 LBL 00
05 ENTER
06 SCI IND Z
07 RND ; round to (z+1) significant digits
08 ISG Z
09 X=Y?
10 GTO 01
11 RDN
12 GTO 00
13 LBL 01
14 RCL Z
15 INT
16 FIX 4 ; optional ;)
.END.
Dieter
Great alternative challenge and greater solution!
Gerson.
I will keep this one. Thanks!
Gerson.
Here are my 2 versions for the WP34s:
Direct version:

3
1/x
RCL[times] L
DEC X
EXPT
+/
RTN
Loop version:

[cmplx]# 1
INC Y
SDR 01
INC X
FP
x=0?
BACK 005
DROP
RTN
Here what I have in my hp 50g, based on Werner's:
<< .6 DUP INV INV  XPON NEG >>
MANT, EXPT— WP 34S has them all. What more is missing? :)
Do you mean SDIG? Aside from solving Dieter's challenge in one step is there any other application for this one? :)
I can't find either of these two in the AUR. Are these really stock functions?
ND1 extends SIZE to Reals.
So #1 can be written
<< pi SIZE 1 >> // machine resolution for Reals
or
<< SIZE(pi)1 >> // if you like algebraic...
and #2 can be written
<< MANT SIZE 1  >> // # significant digits for input Real
Edited: 18 Jan 2012, 12:36 p.m.
For the record: here's one for more recent models featuring a Floor function. Even though it may be called INTG. ;)
It also handles x=0 differently so that the routine quits at one single common endpoint.
M001 LBL M
M002 ENTER
M003 ABS
M004 X<>0?
M005 LOG
M006 INTG
M007 10^x
M008 /
M009 RTN
And finally, here's a completely different approach for the HP41 with XFunctions resp. the 41CX. The goal was: Just display the mantissa, don't disturb the stack and leave everything as it is, including the display format.
The trick here is this: a tendigit mantissa, one decimal marker and one leading space or minussign sum up to twelve alpha characters which will just fit in two registers, thus cutting off a trailing exponent such as "E4".
01 LBL"MANT"
02 STO 01
03 " "
04 X<0?
05 CLA
06 CLX
07 RCLFLAG
08 SCI 9
09 ARCL 01
10 STOFLAG
11 CLX
12 RCL 01
13 ASTO 01
14 ASHF
15 ASTO 02
16 CLA
17 ARCL 01
18 ARCL 02
19 AVIEW
20 END
Yes, users with synthetic programming skills may provide a more elegant solution without any register use. ;)
Dieter
Oliver, we are talking about the WP 34S function set, not RPL.
Oops, sorry, I misread this (it seemed to me in response to Gerson's RPL code). Thanks for the clarification.
Here is a version that works with binary arithmetic, such as Free42/Binary. The results is not exactly an integer since binary does not handle decimal digits as such.
Free42/Binary for Windows returns 15.95.
As it is less obvious than the BCD version, I added some comments.
01 LBL "T"
02 1.5
03 1/X ; this creates a mantissa with 101010...
04 ENTER
05 ENTER
06 1
07 + ; shifts the mantissa 1 bit right, loosing last bit
08 1
09  ; shifts back the mantissa left
10  ; difference is +/1 LSB of mantissa
11 ABS
12 LOG
13 +/
14 END
Nice one! The stringbased ones feel like cheats. I like that this uses purely arithmetic to figure the machine resolution.
Translated to RPL, this returns 15.95, too, on my machine.
Just for completeness: this method fails on a machine with binary arithmetic. Another instance of IEEE754 cleverness. Doublereciprocation of 0.6 yields 0.6.
The 34S doesn't have anything corresponding to the 71's NEIGHBOR function. Well, yet at least :)
Anyone care to do a user code implementation of this one? I think it can be done with a fairly short piece of code.
 Pauli
Edited: 18 Jan 2012, 4:17 p.m.
I'm not quite sure whether the 34s needs a neighbor function, but I think a function that returns 1 ULP of X might be useful. In this case, neighbor(x) simply equals x plus or minus ulp(x).
Assuming standard 16digit precision, ULP(x) simply is 10^(EXPT(x)15):
EXPT
1
5

10^x
Add an ENTER before the first line and a plus or minus at the end, and voilà, you got your neighbor function.
The "15" could be replaced by a special function that returns the current working precision (16 or 34). And if SDL and SDR would accept negative arguments they could be used instead of the final 10^x. Now, if I only knew the English term for "Wink mit dem Zaunpfahl"... ;)
Dieter
Quote:
The "15" could be replaced by a special function that returns the current working precision (16 or 34).
We've got a conditional test to determine which mode the device is in. The ULP function should also work in integer mode returning '1'. Then that are infinities and NaN and subnormals to consider  things are never as easy as they look :(
Quote:
And if SDL and SDR would accept negative arguments they could be used instead of the final 10^x.
This works I think:
ABS
SDR or SDL
Maybe throw in a conditional and a branch to work either way.
Negative arguments are difficult given the way the opcodes are currently structured. Not impossible, just difficult.
Quote:
Now, if I only knew the English term for "Wink mit dem Zaunpfahl"... ;)
Google translate gives "hint hint" which seems appropriate but I'm not sure if the intended meaning is quite the same.
 Pauli
Quote:
Now, if I only knew the English term for "Wink mit dem Zaunpfahl"... ;)
The term you're probably looking for is "broad hint".
For the Germanistically challenged, the literal translation of "Wink mit dem Zaunpfahl" is "hint with a fencepost".
? that has nothing to do with cleverness, but with the different base and smaller ulp.
okay, then this:
.82
ENTER
1/X
1/X

ABS
LOG
CHS
works on 41C, 42S, Free42 Decimal and Binary
Werner
Quote:
can't open the link
Yes, it's once again not working! :(
This crappy '110mb' is the most unreliable webspace server I've ever seen  you should try it again a few times in the next days.
Quote:
? that has nothing to do with cleverness, but with the different base and smaller ulp.
Ok, thanks, I can understand that explanation and confirm that
<< .82 DUP INV INV  ABS LOG NEG >>
returns 15.95... on my machine too.
Sorry for using the word "cleverness". I'm inept at anything related to numerical analysis, and had guessed that maybe there was a design goal to build 1/x in IEEE754 so that 1/(1/x) == x. Upon reflection, that requires infinite (and thereby theoretical) machine resolution.
In gdb:
p 1.0/(1.0/0.6) == 0.6 ==> true
p 1.0/(1.0/0.82) == 0.82 ==> false
Interestingly (to the goof I am), both 0.6 and 0.82 start out as nonexactly representable numbers in IEEE754.
Can you tell if a given number (or its inverse) will pass this test, or not, without trying it out?
In trying this out, I found an error in Free42 Binary (well, rather an 'inexactitude' ;)
In Free42 Binary, do:
0.11 ENTER 11 ENTER 100 / 
it is not zero, and it should be.
The 11 ENTER 100 / is correct (well to the full 52+1 bits), but the 0.11 constant is not, it misses the last bit (in the 52nd position, with the first bit implied).
No doubt this happens because the decimaltobinary conversion just stops at bit 52  truncating the result instead of rounding it, as 11 ENTER 100 / does.
in HEX,
0.11d = 0.1C28F5C28F5C28h
11d/100d = 0.1C28F5C28F5C29h
In Free42 Decimal, on the other hand, no roundtoeven occurs (while it is of course implemented in the real 42S).
Still, I bow in wonder to the amazing achievement of Thomas Okken.
Werner
Here's my first attempt. Sorry, no assemblercompliant output. ;)
01 LBL"ULP"
02 1
03 INTM?
04 SKIP 014
05 DROP
06 SPEC?
07 ERR 01
08 ENTER
09 EXPT
10 3
11 3
12 DBL?
13 SKIP 003
14 DROP
15 1
16 5
17 
18 10^x
19 X<>Y
20 SIGN
21 STOx Y
22 DROP
23 END
The routine returns 1 ULP of x both in integer and real mode. The result is signed, i.e. the output has the same sign as the input. LastX (err... register L) is set correctly. All "special" arguments like NaN or infinity throw a domain error.
Regarding SDL/SDR and the mysterious "Wink mit dem Zaunpfahl": I realized that these shift functions are not very useful here since their arguments cannot exceed 99. So they will not work for very large or small arguments. A good old 10^x will do the trick.
There is no wrapper that saves and restores the stack, and local variables are not used either. I simply have not looked into this yet.
Dieter
Don't worry about the stack mechanics, we can take care of this with some XROM only commands.
There is no need to worry about stack dynamics or last X, we've special commands for dealing with these automatically. Writing proper internal functions in user code is easy now :)
We can extend the range of SDL/SDR arguments to 254 but that still won't help since exponents go much further.
 Pauli
Pauli, indirect arguments may well have much higher limits. It might need a modification but not a major one.
Quote:
can't open the link
Right now it's working again, so take the chance ... :)
For many commands having larger indirect arugments simply doesn't make sense. e.g. FIX, WSIZE, RL, SL, STO, ... The number of commands where a larger limit is sensible is really quite small and I can't think of any where no limit is good.
There are a bit or two left in the argument command table and we could allow an increased limit for indirect but restrict direct arguments to 254. Sounds potentially confusing to me.
 Pauli
The 34S will have the ULP function in the next firmware drop.
 Pauli
Here is a slightly simplified version. It also does not take special care of LastX:
01 LBL"ULP"
02 1
03 INTM?
04 SKIP 013
05 DROP
06 SPEC?
07 ERR 01
08 ENTER
09 EXPT
10 1
11 8
12 DBL?
13 STO+ X
14 
15 3
16 +
17 10^x
18 X<>Y
19 SIGN
20 *
21 END
I don't know if the 34s has a errorignore flag (like the 41's Flag 25). If there is such a flag, the error command in line 07 also has to quit the function. This may be done implicitely (user code functions automatically terminate as soon as an error occurs), or two more lines SPEC? RTN are added.
BTW  I assume the original reason for the discussion on SDL/SDR instead of 10^x was execution speed. I wonder if 10^x really is so much slower that the additional code required for a SDL/SDR solution still makes the whole thing faster. On the other hand I did not look into the Ccode for 10^x. There might be some more or less substantial improvement speedwise if the 10^x code first checks whether x is integer or not. In the first case the result simply hat a mantissa of 1 and an exponent of x  no further mathematics required.
BTW2  in other 34s code samples I noticed constants like #0 or #1. What's the use of these, and in which way do they differ from a simple 0 or 1? Would it make sense to use a "#1" in line 02 here?
Dieter
# 0 differs from 0 in that no interpretation of the input has to take place. Think of # as a special RCL command that addresses a list of constants instead of memory registers.
This doesn't seem to work in integer mode ("bad mode error"). Is this a bug or a feature?
Dieter
Quote:
This doesn't seem to work in integer mode ("bad mode error"). Is this a bug or a feature?
A 'feature' I'm afraid. We have stored the constants in the internal BCD register format only, not as integers. There is no implied conversion (which would be costly in terms of processor cycles, more than interpreting a digit in the code).
Plus only two of the constants are integers, so supporting such conversion doesn't seem like a big benefit.
 Pauli
A few questions:
What is the result for ULP(0) ? 0, the smallest normalised number or the smallest denormalised number? I'll likely switch to the smallest positive denormalised number instead of 0.
Why the sign preservation?
For integer mode I've used a separate routine:
STO L
CLx
INC X
RTN
The function table dispatches separately for integer and real modes so this is free space wise. For real mode, I've decided to code this directly in C  I can play with the exponent field directly and thus avoid the 10^x. Also, handling denomalised numbers properly is something that I think will be difficult for keystroke programs to deal with.
 Pauli
If we want to preserve sign and return 0 for 0 the function is identical to SIGN in integer mode. With Pauli's implementation of ULP (any mode), multiplication with SIGN will return the result of Dieter's routine.
Quote:
A few questions:
What is the result for ULP(0) ? 0, the smallest normalised number or the smallest denormalised number?
For the HP71B, the NEIGHBOR function returns the smallest positive denormalized number if the second argument is positive, the smallest negative denormalized number if the second argument is negative, and 0 if the second argument is 0.
If you haven't set the flag to allow for denormalized results, you'll get an Underflow warning and a default 0 result instead.
Best regards from V.
I think MINVAL & MAXVAL (mode dependent) may be useful information. Denormalized numbers are always allowed in WP 34S, only NaNs and infinities need a system flag to be set to not cause an error.
Quote:
I'm not quite sure whether the 34s needs a neighbor function, but I think a function that returns 1 ULP of X might be useful. In this case, neighbor(x) simply equals x plus or minus ulp(x).
Actually it doesn't work this way :(
ULP(10) = 1e14. 10  1e14 = 9.99999999999999.
The correct neighbour below 10 is 1e15 below = 9.999999999999999
i.e. 10  ULP(10) loses a digit.
 Pauli
First, it's easier to understand if these numbers with lots of zeros or nines are written with grouped digits. ;)
Yes, exact powers of 10 are a special case. The next higher number here actually is 10 + 1E14:
10 + 1E14 = 10,00 0000 0000 0001
Which is what the proposed routine returns: ULP(10) = 1E14.
On the other hand, the next lower value is not 1E14 but just 1E15 below 10:
10  1E15 = 9,9999 9999 9999 9999
So there are actually two different ways to handle this.
This could be solved by two separate ULP functions, one for each case (up/down). Or we do it the same way as the 71B and implement a neighbor function with two arguments. Or... ;)
Dieter
I was copy/pasting the numbers from the console emulator that doesn't include digit separators or grouping :)
ULP is good the way it is. It returns 1 ULP. Splitting this for ULP upwards and ULP downwards seems wrong to me. Having a neighbour function like the 71B is one option, having a pair of functions that return the next value larger or smaller than their argument is another (which I think I prefer but not by much).
Anyway, I've got ULP working properly for subnormals. For subnormals (and zero) it returns the smallest subnormal.
 Pauli
With regard to the neighbo(u)r function:
Quote:
...a pair of functions that return the next value larger or smaller than their argument is another (which I think I prefer but not by much).
Then throw in my preference for the latter solution. In addition to a (signed or unsigned) ULP function. Two separate neighbor functions are fine if they return two values that are not larger resp. smaller, but closer to zero resp. farther out. This difference is important for negative arguments. That's also the main reason why my original ULP routine returns a signed result: adding resp. subtracting this to the original X consistently returns the one or the other.
Dieter
Quote:
No doubt this happens because the decimaltobinary conversion just stops at bit 52  truncating the result instead of rounding it, as 11 ENTER 100 / does.
I think you're right. I haven't touched that code in years and it's a bit light on comments, so I'll have to stare at it a bit more to make sure. The code is *supposed* to keep going through successively lower powers of 2 until they become too small to affect the mantissa, so it should give the correct result, but since it doesn't, I'll have to dust off the debugger and find out why...
 Thomas
I think it would thus be better to have four functions which give the next representable value in the direction:
towards 0
away from 0
towards negative infinity
towards positive infinity
Maybe the 71b's neighbour or C nextafter function is what we want here after all. Then the user can specify exactly which direction to go.
 Pauli
For this, I can use a workaround  just enter the constants in a different way (like 11 ENTER 100 /), but for the Decimal version not performing a roundtoeven is much more fundamental to me, and there's no way round:
3 1/X 2 /
exact: 0.1666...6665, where the trailing 5 will be rounded.
42S : 0.166666666666 (roundtoeven)
Free42 Dec: 0.16666 66666 66666 66666 66667 (round up)
Or, for addition:
42S : 1 ENTER 5e12 + > 1.00000000000
Free42 Dec: 1 ENTER 5e25 + > 1.00000 00000 00000 00000 0001
(all tested on iPhone Free42 version)
Werner
Edited: 22 Jan 2012, 3:55 a.m.