Talk:Arithmetic-geometric mean: Difference between revisions
Content added Content deleted
m (→rapidity: added a new talk section.) |
m (→rapidity: added a new talk section, fixed/(supplied) the previous LANG ending HTML tag.) |
||
| Line 47: | Line 47: | ||
) |
) |
||
const ε = 1e-14 |
const ε = 1e-14</lang> |
||
==rapidity |
==rapidity== |
||
From this Rosetta Code task's prologue: |
|||
func agm(a, g complex128) complex128 { |
|||
for cmplx.Abs(a-g) > cmplx.Abs(a)*ε { |
|||
a, g = (a+g)*.5, cmplx.Rect(math.Sqrt(cmplx.Abs(a)*cmplx.Abs(g)), |
|||
(cmplx.Phase(a)+cmplx.Phase(g))*.5) |
|||
} |
|||
return a |
|||
} |
|||
:: <big><big>Since the limit of <math>a_n-g_n</math> tends (rapidly) to zero with iterations, this is an efficient method.</big></big> |
|||
func main() { |
|||
fmt.Println(agm(1, 1/math.Sqrt2)) |
|||
}</lang> |
|||
<br>With this in mind, I modified the REXX's entry (and suppressed the normal output being displayed), and |
|||
End of copy --[[User:Sonia|Sonia]] ([[User talk:Sonia|talk]]) 19:16, 17 June 2013 (UTC) |
|||
<br>added a display of the iteration count along with the number of decimal digits being used. |
|||
The modified REXX program is: |
|||
<lang rexx>/*REXX program calculates the AGM (arithmetic─geometric mean) of two numbers.*/ |
|||
parse arg a b digs . /*obtain optional numbers from the C.L.*/ |
|||
if digs=='' | digs==',' then digs=100 /*No DIGS specified? Then use default.*/ |
|||
numeric digits digs /*REXX will use lots of decimal digits.*/ |
|||
if a=='' | a==',' then a=1 /*No A specified? Then use default.*/ |
|||
if b=='' | b==',' then b=1/sqrt(2) /*No B specified? " " " */ |
|||
call AGM a,b /*invoke AGM & don't show A,B,result.*/ |
|||
exit /*stick a fork in it, we're all done. */ |
|||
/*────────────────────────────────────────────────────────────────────────────*/ |
|||
agm: procedure: parse arg x,y; if x=y then return x /*equality case?*/ |
|||
if y=0 then return 0 /*is Y zero? */ |
|||
if x=0 then return y/2 /* " X " */ |
|||
d=digits(); numeric digits d+5 /*add 5 more digs to ensure convergence*/ |
|||
tiny='1e-' || (digits()-1); /*construct a pretty tiny REXX number. */ |
|||
ox=x+1 |
|||
do #=1 while ox\=x & abs(ox)>tiny; ox=x; oy=y |
|||
x=(ox+oy)/2; y=sqrt(ox*oy) |
|||
end /*while ··· */ |
|||
numeric digits d /*restore numeric digits to original.*/ |
|||
/*this is the only output displayed ►─┐*/ |
|||
say 'digits='right(d,7)", iterations=" right(#,3) /* ◄───────────────┘*/ |
|||
return x/1 /*normalize X to the new digits. */ |
|||
/*────────────────────────────────────────────────────────────────────────────*/ |
|||
sqrt: procedure; parse arg x; if x=0 then return 0; d=digits(); i=; m. =9 |
|||
numeric digits 9; if x<0 then do; x=-x; i='i'; end; numeric form; h=d+2 |
|||
parse value format(x,2,1,,0) 'E0' with g 'E' _ .; g=g*.5'E'_%2 |
|||
do j=0 while h>9; m.j=h; h=h%2+1; end /*j*/ |
|||
do k=j+5 to 0 by -1; numeric digits m.k; g=.5*(g+x/g); end /*k*/ |
|||
numeric digits d; return (g/1)i</lang> |
|||
'''outputs''' with a specified number of digits, and with the various outputs being combined/reduced/edited: |
|||
<pre> |
|||
digits= 100, iterations= 9 |
|||
digits= 200, iterations= 10 |
|||
digits= 400, iterations= 11 |
|||
digits= 800, iterations= 12 |
|||
digits= 1600, iterations= 13 |
|||
digits= 3200, iterations= 14 |
|||
digits= 6400, iterations= 15 |
|||
digits= 12800, iterations= 16 |
|||
digits= 25600, iterations= 17 |
|||
digits= 51200, iterations= 18 |
|||
</pre> |
|||
Needless to say, the CPU time increased everytime the digits were doubled. |
|||