The Museum of HP Calculators

HP Forum Archive 21

 Need help understanding math....Message #1 Posted by cyrille de Brébisson on 9 Dec 2013, 1:40 p.m. Hello, I am trying to follow up on a question that I have had for a long time on the HP48 code for the complex acos and asin functions... The comment says to test against 1/sqrt(2), but the code tests against 0.7! It looks to me that the 2 cases are there to improve precision in upper and lower octants, and that having the midpoint a little bit off is of little to no consequences, but I wanted to see what the math guru in the forum thought about it. cyrille * Q := y^2/(sqrt((|x|+1)^2 + y^2) + (|x|+1)) * R := sqrt((|x|-1)^2 + y^2) +||x|-1| * S := y^2/R if R<>0, 0 otherwise * M := Q+R if |x|>=1, Q+S otherwise * P := Q+S if |x|>=1, Q+R otherwise * B := 2*x/(M+2) * C := sqrt((P/(M+2))*(2-(P/(M+2)))) * sg(y,x) := sgn(y) if y<>0, -sgn(x) otherwise * IM := sg(y,x)*lnp1((M/2) + sqrt((M/2)*((M/2)+2))) (sign replacement) * * { arccos(B) |B| <= (1/sqrt(2)) * RE1 := { arcsin(C) |B| > (1/sqrt(2)) and B >= 0 * { pi - arcsin(C) |B| > (1/sqrt(2)) and B < 0 * * RE2 := { arcsin(B) |B| <= (1/sqrt(2)) * { sgn(B)*arccos(C) |B| > (1/sqrt(2)) (sign replacement) * :: C%>%% 2DUP DUP %%*SWAP %%ABS DUP %%1+ DUPDUP %%* 4PICK %%+ %%SQRT %%+ 3PICKSWAP %%/ OVER %%1 %%- %%ABS DUPDUP %%* 5PICK %%+ %%SQRT %%+ 4ROLLOVER DUP %%0<> ITE %%/ :: 2DROP %%0 ; UNROTOVER %%+ UNROT %%+ ROT %%1 %%< ?SKIPSWAP DUP %%2 %%+ ROTOVER %%/ %%2 OVER %%- %%* %%SQRT 5PICK DUP %%+ ROT %%/ ROT %%2 %%/ DUP %%2 %%+ OVER %%* %%SQRT %%+ %%LNP1 5ROLL 5ROLL DUP %%0= ITE :: DROP %%0< ; :: SWAPDROP %%0>= ; ?SKIP %%CHS UNROTDUP %%ABS %%.7 %%<= case :: SWAPDROPDUP %%ACOSRAD SWAP %%ASINRAD ; SWAPDUP %%ASINRAD SWAP %%ACOSRAD %%ABS ROT %%0>= ?SEMI %%CHS %%PI ROT %%- SWAP ; Edited: 9 Dec 2013, 1:46 p.m.

 Re: Need help understanding math....Message #2 Posted by Han on 9 Dec 2013, 2:04 p.m.,in response to message #1 by cyrille de Brébisson Edited for legibility: Quote: Hello, I am trying to follow up on a question that I have had for a long time on the HP48 code for the complex acos and asin functions... The comment says to test against 1/sqrt(2), but the code tests against 0.7! It looks to me that the 2 cases are there to improve precision in upper and lower octants, and that having the midpoint a little bit off is of little to no consequences, but I wanted to see what the math guru in the forum thought about it. cyrille ```* Q := y^2/(sqrt((|x|+1)^2 + y^2) + (|x|+1)) * R := sqrt((|x|-1)^2 + y^2) +||x|-1| * S := y^2/R if R<>0, 0 otherwise * M := Q+R if |x|>=1, Q+S otherwise * P := Q+S if |x|>=1, Q+R otherwise * B := 2*x/(M+2) * C := sqrt((P/(M+2))*(2-(P/(M+2)))) * sg(y,x) := sgn(y) if y<>0, -sgn(x) otherwise * IM := sg(y,x)*lnp1((M/2) + sqrt((M/2)*((M/2)+2))) (sign replacement) * * { arccos(B) |B| <= (1/sqrt(2)) * RE1 := { arcsin(C) |B| > (1/sqrt(2)) and B >= 0 * { pi - arcsin(C) |B| > (1/sqrt(2)) and B < 0 * * RE2 := { arcsin(B) |B| <= (1/sqrt(2)) * { sgn(B)*arccos(C) |B| > (1/sqrt(2)) (sign replacement) * :: C%>%% 2DUP DUP %%*SWAP %%ABS DUP %%1+ DUPDUP %%* 4PICK %%+ %%SQRT %%+ 3PICKSWAP %%/ OVER %%1 %%- %%ABS DUPDUP %%* 5PICK %%+ %%SQRT %%+ 4ROLLOVER DUP %%0<> ITE %%/ :: 2DROP %%0 ; UNROTOVER %%+ UNROT %%+ ROT %%1 %%< ?SKIPSWAP DUP %%2 %%+ ROTOVER %%/ %%2 OVER %%- %%* %%SQRT 5PICK DUP %%+ ROT %%/ ROT %%2 %%/ DUP %%2 %%+ OVER %%* %%SQRT %%+ %%LNP1 5ROLL 5ROLL DUP %%0= ITE :: DROP %%0< ; :: SWAPDROP %%0>= ; ?SKIP %%CHS UNROTDUP %%ABS %%.7 %%<= case :: SWAPDROPDUP %%ACOSRAD SWAP %%ASINRAD ; SWAPDUP %%ASINRAD SWAP %%ACOSRAD %%ABS ROT %%0>= ?SEMI %%CHS %%PI ROT %%- SWAP ; ```

 Re: Need help understanding math....Message #3 Posted by cyrille de Brébisson on 9 Dec 2013, 3:17 p.m.,in response to message #2 by Han thanks, I tried 4 times but could not get it to take my CR properly! cyrille

 Re: Need help understanding math....Message #4 Posted by Walter B on 9 Dec 2013, 3:29 p.m.,in response to message #3 by cyrille de Brébisson Once read that rats learn it with four times ;-) SCNR

 Re: Need help understanding math....Message #5 Posted by Jeff O. on 9 Dec 2013, 3:42 p.m.,in response to message #3 by cyrille de Brébisson You just need to enclose your formatted text between [pre] and [/pre] tags.

 Re: Need help understanding math....Message #6 Posted by David Hayden on 12 Dec 2013, 5:28 p.m.,in response to message #2 by Han I don't understand the comments in this code. What are Q, R, S, M, P, B, and C? Is RE1 the real part of stack level 1? RE2 is real part of stack level 2? IM is imaginary part of... both? Just curious. Dave

 Re: Need help understanding math....Message #7 Posted by Thomas Klemm on 16 Dec 2013, 6:59 p.m.,in response to message #6 by David Hayden HTH Thomas

 Re: Need help understanding math....Message #8 Posted by Paul Dale on 9 Dec 2013, 11:44 p.m.,in response to message #1 by cyrille de Brébisson I've not put any effort into investigating this and am really just guessing, however I suspect that your intuition here is correct. That area of the sine and cosine functions is pretty smooth and well behaved -- essentially linear over the reals and I don't think it goes weird as one moves away from the real line. I'd guess it was a case of 0.7 being an existing constant that was near enough. Depending on how the comparison is coded, it might be very slightly faster too but unlikely enough to notice. - Pauli

 Re: Need help understanding math....Message #9 Posted by Raymond Del Tondo on 12 Dec 2013, 5:42 p.m.,in response to message #8 by Paul Dale Quote:I'd guess it was a case of 0.7 being an existing constant that was near enough. Depending on how the comparison is coded, it might be very slightly faster too but unlikely enough to notice.Actually the constant is defined in that code slice, and not a call to a constant defined elsewhere, thus saving one pointer resolution time cycle. I'm confident the common code for complex ASIN and ACOS dates back to the Paladin, so could simply be a LISARIF (Bill Wickes)

 Re: Need help understanding math....Message #10 Posted by Didier Lachieze on 13 Dec 2013, 2:23 a.m.,in response to message #9 by Raymond Del Tondo Quote:I'm confident the common code for complex ASIN and ACOS dates back to the Paladin, so could simply be a LISARIF (Bill Wickes) Some decoding for people who - like me - would need it: Paladin: internal hp code name for the HP-28C LISARIF: "Life is Short and ROM is Full" a famous quote from Bill Wickes (see the original source)

Go back to the main exhibit hall