New guy and programming problem

[Jlouis]

(08-27-2017 07:52 PM)pier4r Wrote:  

(08-27-2017 06:11 PM)Dieter Wrote:  Everyone knows that 1+1=2 but your calculator does it anyway. ;-)

Dieter

Until you buy a gaxio.

I enjoy this review very much. It's early in the morning and I have already done good laughs.

Cheers

[Logan]

So for completeness' sake, I'll add the response I left to you when you asked this question on my 42s programming video.

Sure, so what's happening is that if you enter a new value for R1 or R2 the stack is lifting and the old value is now in the Y portion of the stack, so ultimately your program is calculating the old R2 + 1/(new R2) and then 1/x.

A couple of solutions: the most straightforward is probably to store the intermediate value, so after the first 1/X, STO 00 and then after the second 1/X, RCL 00.

Or, perhaps a better practice is to create the variables first, then perform the operations:

Code:

LABEL "PAR2"
INPUT R1
INPUT R2
1/X
RCL R1
1/X
+
1/X
END

Lastly, if you wanted to be really streamlined, you could do away with INPUT altogether and just enter your values into the stack before executing the program

Code:

LABEL "PAR2"
1/X
X <> Y
1/X
+
1/X
END

(Run by entering R1 ENTER R2 XEQ PAR2)

[Skirmitch]

(08-28-2017 02:20 PM)Logan Wrote:  So for completeness' sake, I'll add the response I left to you when you asked this question on my 42s programming video.

Sure, so what's happening is that if you enter a new value for R1 or R2 the stack is lifting and the old value is now in the Y portion of the stack, so ultimately your program is calculating the old R2 + 1/(new R2) and then 1/x.

A couple of solutions: the most straightforward is probably to store the intermediate value, so after the first 1/X, STO 00 and then after the second 1/X, RCL 00.

Or, perhaps a better practice is to create the variables first, then perform the operations:

Code:

LABEL "PAR2"
INPUT R1
INPUT R2
1/X
RCL R1
1/X
+
1/X
END

Lastly, if you wanted to be really streamlined, you could do away with INPUT altogether and just enter your values into the stack before executing the program

Code:

LABEL "PAR2"
1/X
X <> Y
1/X
+
1/X
END

(Run by entering R1 ENTER R2 XEQ PAR2)

Yes, I did read your reply, also lemme use the instance to thank you so much for your videos, I've watched every single one of them (calculators related) and its partially your fault that I'm so inclined to the 42S hahaha. Greetings

[Dieter]

(08-27-2017 11:12 PM)Paul Dale Wrote:  Yes to your other suggestions for improving accuracy, I got in what I could -- limited flash space and running out of energy justifying everything.

Life is short and RAM is full. ;-)

(08-27-2017 11:12 PM)Paul Dale Wrote:  Apart from 1 - sin(x) which doesn't have accuracy concerns that I can see

Well, for 90±5E–7° even on the 34s sin(x) rounds to 1 (in SP mode).
On a 10-digit machine this already happens for 89,9995°...90,0005°.
Even for 89°...91° three significant digits of 1 – sin(x) are lost.

(08-27-2017 11:12 PM)Paul Dale Wrote:  ...and √(x²+y²) which is present already as complex ABS.

Ah, yes, that's something I forgot. I even use complex ABS for this purpose here and there, both on the 34s and the 35s.

Dieter

[Logan]

Ha! I'll accept that blame!

[Paul Dale]

(08-28-2017 06:13 PM)Dieter Wrote:  Well, for 90±5E–7° even on the 34s sin(x) rounds to 1 (in SP mode).
On a 10-digit machine this already happens for 89,9995°...90,0005°.
Even for 89°...91° three significant digits of 1 – sin(x) are lost.

Good point, I was wrong. For some reason I thought the flattest gradient was 1 not 0.


Pauli

[Thomas Klemm]

(08-28-2017 06:13 PM)Dieter Wrote:  On a 10-digit machine this already happens for 89,9995°...90,0005°.

We can use this difference-to-product formula to calculate e.g. \(1-\sin(89.9995°)=\sin(90°)-\sin(89.9995°)\):

\(\sin \theta - \sin \varphi =2\sin \left({\frac {\theta - \varphi }{2}}\right)\cos \left({\frac {\theta + \varphi }{2}}\right)\)

And then since \(\theta = \frac{\pi}{2}\) we can set \(\varepsilon = \frac {\theta - \varphi }{2}\) and use \(\cos \left({\tfrac {\pi }{2}}-\varepsilon \right)=\sin \varepsilon\) and end up with:

\(1 - \sin \varphi =2\sin^2 \varepsilon\)

Thus we set \(\theta=90\) and \(\varphi=89.9995\) and get:

90
ENTER
89.9995
-
2
÷
SIN
x²
2
×

3.807717748e-11

(Calculated using an HP-41.)

Cheers
Thomas

[Albert Chan]

(09-01-2018 12:23 PM)Thomas Klemm Wrote:  

(08-28-2017 06:13 PM)Dieter Wrote:  On a 10-digit machine this already happens for 89,9995°...90,0005°.

We can use this difference-to-product formula to calculate e.g. \(1-\sin(89.9995°)=\sin(90°)-\sin(89.9995°)\) ...

Another way, rearrange 1 - sin(x), so more accurate for x close to right angle.

1 - sin(x)
= (1 - sin(x)) * (1 + sin(x)) / (1 + sin(x))
= cos(x)^2 / (1 + sin(x))

1 - sin(89.9995°) ~ cos(89.9995°)^2 / 2 ~ 3.807717747e-11

[sasa]

(09-01-2018 12:23 PM)Thomas Klemm Wrote:  We can use this difference-to-product formula to calculate...
...
3.807717748e-11

Even result is good enough, there is still too much trouble to get accurate result to the last significant digit with arguments near 90 degree.

The first problem is what is each calculator's internal mantissa length, the second exact representation for PI in internal format (internal conversion from degree to radian) and proper round off mechanism. And as well the used format (decimal or binary floating point) and proper algorithms to calculate trig functions.

The Wolfram Alpha returns:
3.8077177473096929677620869980236729697355177044497495... × 10^-11

[Thomas Klemm]

(09-01-2018 03:18 PM)sasa Wrote:  Even result is good enough, there is still too much trouble to get accurate result to the last significant digit with arguments near 90 degree.

For sin(0.00025°) WolframAlpha gives:

4.3633231299719786898764073194250596667396107572873548... × 10^-6

The HP-41C returns:

4.363323130e-06

This result is correctly rounded to the last digit.

Then we square the rounded number and get:

1.90385887367929969 × 10^-11

The HP-41C returns:

1.903858874e-11

Again the result is rounded correctly.

Last thing to do is to multiply the number by 2:

3.807717748 × 10^-11

Which is exactly what the HP-41C returns:

3.807717748e-11


As we can see the HP-41C rounds these intermediate results correctly.

The problem with the cancelation when calculating \(1-sin(x)\) is not related to the accuracy of calculating \(sin(x)\) for \(x\) close to \(90°\).

---

Similar to Albert's trick we can multiply \(\sqrt{1+x^2}-1\) by \(\frac{\sqrt{1+x^2}+1}{\sqrt{1+x^2}+1}=1\) and get:

\(\begin{align*}
\sqrt{1+x^2}-1 &= \sqrt{1+x^2}-1 \cdot \frac{\sqrt{1+x^2}+1}{\sqrt{1+x^2}+1} \\
&= \frac{(\sqrt{1+x^2}-1)(\sqrt{1+x^2}+1)}{\sqrt{1+x^2}+1} \\
&= \frac{1+x^2-1}{\sqrt{1+x^2}+1} \\
&= \frac{x^2}{\sqrt{1+x^2}+1}
\end{align*}\)

And again cancelation for small \(x\) is avoided.

[Albert Chan]

(09-01-2018 03:18 PM)sasa Wrote:  Even result is good enough, there is still too much trouble to get accurate result
to the last significant digit with arguments near 90 degree ...

The Wolfram Alpha returns:
3.8077177473096929677620869980236729697355177044497495... × 10^-11

You might have the impression that Wolfram answer is easily calculated. It is *not*

You don't feel the difficulty, because modern computers are just too fast.
Or, in the old days, where you actually get *charge* for computing time

If you don't believe me, try to reproduce Wolfram answer by hand.

[sasa]

(09-01-2018 05:20 PM)Albert Chan Wrote:  You might have the impression that Wolfram answer is easily calculated. It is *not*

Yes, I have had to go myself through all these headache with my own library for MCUs...

Simply, there is too many things to balance and decide proper path, if resources are quite limited.

[Vtile]

(08-27-2017 05:12 PM)Gerson W. Barbosa Wrote:  Is there another calculator that offers this as a built-in function?

Gerson.

Yes and it is even a standard factory build: Casio FX-61F.

What comes to the OPs original question.. Dieters two last program snippets are the way it should be done with the RPN or RPL machine (atleast the 2ns last should also work without issues with impedances).
Then there is some ways to do it "better" to avoid some insignificant rounding errors which are irrelevant most all practical levels.

[Albert Chan]

Third way to calculate accurate 1 - sin(x), is to remove sine function, replaced with cosine.
We can then pick whatever double angle formula is suitable:

cos(2x) = cos(x)² - sin(x)² = 2 cos(x)² - 1 = 1 - 2 sin(x)²

1 - sin(89.9995°)
= 1 - cos(0.0005°)
= 1 - (1 - 2 sin(0.00025°)²)
= 2 sin(0.00025°)²
= 3.807717747e-11

Note: above same as Klemm difference-to-product formula (post 27), but probably easier to remember.

BTW, you can build multiple angles cosine formula with cosines ... Neat

cos(n x) = 2 cos(x) cos((n-1) x) - cos((n-2) x)

[Thomas Klemm]

Fun fact: These trigonometric identities were used for multiplication and division before logarithms were known.

Cf. Prosthaphaeresis

Quote:Notice how similar the above algorithm is to the process for multiplying using logarithms, which follows these steps: scale down, take logarithms, add, take inverse logarithm, scale up. It's no surprise that the originators of logarithms had used prosthaphaeresis. Indeed the two are closely related mathematically. In modern terms, prosthaphaeresis can be viewed as relying on the logarithm of complex numbers, in particular on Euler's formula:

\(e^{ix}=\cos x+i\sin x\)

I find it interesting that this technique was developed using complex numbers under the hood while they weren't discovered yet.

Cheers
Thomas

[Albert Chan]

I never knew trig function were used as logarithm !
Example: arccos/cos vs ln/exp

0.123 * 0.456
= cos(82.9347°) * cos(62.8707°)
= 1/2 (cos(82.9347° + 62.8707°) + cos(82.9347° - 62.8707°))
= 1/2 (cos(145.8054°) + cos(20.0640°))
= 0.05608822 ~ 0.056088

0.123 * 0.456
= exp(ln(0.123) + ln(0.456))
= exp(-2.09557 + -0.78526)
= exp(-2.88083)
= 0.05608819 ~ 0.056088

But, how can trig function used for division ?

[Thomas Klemm]

(09-02-2018 09:21 PM)Albert Chan Wrote:  But, how can trig function used for division ?

From the linked Wikipedia article:

Quote:To divide, we exploit the definition of the secant as the reciprocal of the cosine.

[Vtile]

(09-02-2018 09:33 PM)Thomas Klemm Wrote:  

(09-02-2018 09:21 PM)Albert Chan Wrote:  But, how can trig function used for division ?

From the linked Wikipedia article:

Quote:To divide, we exploit the definition of the secant as the reciprocal of the cosine.

Reciprocals <3

[Albert Chan]

(09-02-2018 09:40 PM)Vtile Wrote:  Reciprocals <3

Reciprocals have to be below 1, to satisfy arccos domain:

Example: 0.123 / 0.456 = 10 * 0.123 * (1/4.56)

arccos(1/4.56) = arcsec(4.56) = 77.3322°

0.123 / 0.456
= 10 cos(82.9347°) cos(77.3322°)
= 10/2 * (cos(160.2669°) + cos(5.6025°))
~ 0.26974