03-20-2022, 07:28 PM
(03-20-2022 07:39 AM)Ángel Martin Wrote: [ -> ]Sorry about that, here's the missing code in the [INTZL] routine:
Hi Ángel,
Thanks for providing - works like a charm now!
KR
Frido
(03-20-2022 07:39 AM)Ángel Martin Wrote: [ -> ]Sorry about that, here's the missing code in the [INTZL] routine:
Quote:Let \(I_k\) and \(C_k\) denote the areas of regular \(k\)-gons inscribed in and circumscribed around a given circle.
Then for all \(n\), \(I_{2n}\) is the geometric mean of \(I_n\) and \(C_n\), and \(C_2n\) is the harmonic mean of \(I_2n\) and \(C_n\); that is,
\(
\begin{align}
I_{2n}=\sqrt{I_{n}C_{n}}
\end{align}
\)
and
\(
\begin{align}
C_{2n}=\frac{2}{\frac{1}{I_{2n}}+\frac{1}{C_{n}}}.
\end{align}
\)
001 RCL× Y
002 √
003 ||
004 STO+ x
005 RCL L
006 END
(03-21-2022 07:24 AM)Thomas Klemm Wrote: [ -> ]So I thought I give it a try using the || operator of the WP-34S:
(03-21-2022 04:03 PM)Frido Bohn Wrote: [ -> ]How would the || operator translate into e.g., the HP41?
1/X
X<>Y
1/X
+
1/X
RCL× Y
X<>Y
RCL+ L
/
00 { 16-Byte Prgm }
01 RCL× ST Y
02 SQRT
03 X<>Y
04 RCL× ST Y
05 LASTX
06 RCL+ ST Z
07 ÷
08 STO+ ST X
09 X<>Y
10 END
(03-21-2022 05:27 PM)Thomas Klemm Wrote: [ -> ]But I have no idea how the WP-34S calculates the result.
/**************************************************************************/
/* The real parallel operator.
*
* par(x, y) = (x . y) / (x + y) = 1 / ( 1/x + 1/y )
*/
XLBL"PARL"
xIN DYADIC
x=0?
xOUT xOUT_NORMAL
ENTER[^]
RCL+ Z
x=0?
JMP ret_neginf
/
[times]
xOUT xOUT_NORMAL
00 { 6-Byte Prgm }
01 ENTER
02 RCL+ ST Z
03 ÷
04 ×
05 END
(03-21-2022 07:24 AM)Thomas Klemm Wrote: [ -> ]We end up with the same sequence as when using Vieta's formula.
10 N=6 @ A=SIN(PI/N*2)*N/2 @ B=TAN(PI/N)*N
20 A=SQRT(A*B) @ B=B*A/((B+A)*.5) @ N=N+N
30 DISP N,A,B,(A+2*B)/3
40 IF N<100 THEN 20
>RUN
12 3 3.21539030918 3.14359353945
24 3.10582854123 3.15965994211 3.14171614182
48 3.13262861329 3.14608621514 3.14160034786
96 3.13935020305 3.14271459965 3.14159313412
192 3.14103195089 3.14187304998 3.14159268362
(03-21-2022 10:45 PM)Albert Chan Wrote: [ -> ]Code:
>RUN
12 3 3.21539030918 3.14359353945
24 3.10582854123 3.15965994211 3.14171614182
48 3.13262861329 3.14608621514 3.14160034786
96 3.13935020305 3.14271459965 3.14159313412
a b
48 3.13935020304686722 3.14608621513143498
c d
96 3.14103195089050965 3.14271459964536831
pi ~ (640*c*(28*c-139*a)/(c-4*a)+4*b*(3737+14256*b/(d-4*b))+1056*d-2568*a-6453*(a^2*b)^(1/3))/11655
pi ~ 3.141592653589793(44)
(03-24-2022 01:36 AM)Gerson W. Barbosa Wrote: [ -> ]Code:
a b
48 3.13935020304686722 3.14608621513143498
c d
96 3.14103195089050965 3.14271459964536831
pi ~ (640*c*(28*c-139*a)/(c-4*a)+4*b*(3737+14256*b/(d-4*b))+1056*d-2568*a-6453*(a^2*b)^(1/3))/11655
pi ~ 3.141592653589793(44)
(03-26-2022 03:59 PM)Albert Chan Wrote: [ -> ](03-24-2022 01:36 AM)Gerson W. Barbosa Wrote: [ -> ]Code:
a b
48 3.13935020304686722 3.14608621513143498
c d
96 3.14103195089050965 3.14271459964536831
pi ~ (640*c*(28*c-139*a)/(c-4*a)+4*b*(3737+14256*b/(d-4*b))+1056*d-2568*a-6453*(a^2*b)^(1/3))/11655
pi ~ 3.141592653589793(44)
How is formula derived ?
(03-26-2022 03:59 PM)Albert Chan Wrote: [ -> ]r = (b-a)/(b+b) = (1-a/b)/2 = (1-cos(x))/2 = sin(x/2)^2
x/2 = asin(sqrt(r))
...
A = sin(x)/x * pi
pi = A * (2*x/2) / (2*sin(x/2)*cos(x/2))
= A/sqrt(1-r) * asin(sqrt(r)) / sqrt(r)
= A/sqrt(1-r) * (1 + r/6 * (1 + r*3^2/(4*5) * (1 + r*5^2/(6*7) * (1 + r*7^2/(8*9) + ...
(05-31-2021 07:30 PM)Albert Chan Wrote: [ -> ]\(\displaystyle\arcsin(x) = 2\arcsin\left( {x \over \sqrt{2\sqrt{1-x^2}+2}} \right)\)
function asinq(x) -- = asin(sqrt(x))
if x > 4e-4 then return 2 * asinq(0.5*x/(sqrt(1-x)+1)) end
return sqrt(x) * (1+x/6*(1+x*9/20/(1-x*25/42)))
end
function RC(x,y, verbose) -- RC(1-x,1) = asin(sqrt(x))/sqrt(x)
local k = y + 2*sqrt(x*y)
if k==x+2*y then return sqrt(3/k) end
if verbose then print(sqrt(3/k)) end
return RC((x+k)*0.25, (y+k)*0.25, verbose)
end
function C(x) -- = asin(sqrt(x)) / sqrt(x)
if x < 4e-4 then return 1+x/6*(1+x*9/20/(1-x*25/42)) end
local y = 0.5/(sqrt(1-x)+1)
y = 2*sqrt(y) * C(x*y)
print('C( '..x..' ) = ',y)
return y
end
lua> C(1)
C( 0.0006022718974138037 ) = 1.0001004058641836
C( 0.002407636663901557 ) = 1.0004017081549652
C( 0.009607359798384776 ) = 1.0016081890839748
C( 0.038060233744356624 ) = 1.0064545427995637
C( 0.14644660940672624 ) = 1.0261721529770307
C( 0.5 ) = 1.1107207345395913
C( 1 ) = 1.5707963267948963
(03-19-2022 11:18 AM)Ángel Martin Wrote: [ -> ]Heresy or new truth?
Here's an interesting reading for y'all:
https://www.theverge.com/tldr/2018/3/14/...iday-truth
(03-14-2022 11:06 PM)ttw Wrote: [ -> ]One can often apply various sequence transformations to speed up these series. A couple of new methods are given here.Interesting. When I see PI, I am remembering of page 51 of the "exact" PI formula
https://arxiv.org/pdf/1702.07199.pdf
https://people.mpim-bonn.mpg.de/zagier/f...lltext.pdf
In most record-setting attempts, the desired result is a fraction with really big integers, even bigger than 355/113. An (to me, anyway) amusing feature is that most of the time is taken by the final long division. Third order, first order, fourth order, etc., methods all take about the same amount of time as the time for arithmetic grows pretty fast.