

 This page uses content from Wikipedia. The original article was at Arithmetic-geometric mean. The list of authors can be seen in the page history. As with Rosetta Code, the text of Wikipedia is available under the GNU FDL. (See links for details on variance)



Task

Write a function to compute the arithmetic-geometric mean of two numbers.



The arithmetic-geometric mean of two numbers can be (usefully) denoted as 




a
g
m

(
a
,
g
)


{\displaystyle \mathrm {agm} (a,g)}

, and is equal to the limit of the sequence:

 




a

0


=
a
;


g

0


=
g


{\displaystyle a_{0}=a;\qquad g_{0}=g}


 




a

n
+
1


=



1
2



(

a

n


+

g

n


)
;


g

n
+
1


=



a

n



g

n




.


{\displaystyle a_{n+1}={\tfrac {1}{2}}(a_{n}+g_{n});\quad g_{n+1}={\sqrt {a_{n}g_{n}}}.}


Since the limit of 




a

n


−

g

n




{\displaystyle a_{n}-g_{n}}

 tends (rapidly) to zero with iterations, this is an efficient method.

Demonstrate the function by calculating:






a
g
m

(
1
,
1

/



2


)


{\displaystyle \mathrm {agm} (1,1/{\sqrt {2}})}




Also see

   mathworld.wolfram.com/Arithmetic-Geometric Mean


-- ARITHMETIC GEOMETRIC MEAN -------------------------------------------------
 
property tolerance : 1.0E-5
 
-- agm :: Num a => a -> a -> a
on agm(a, g)
    script withinTolerance
        on |λ|(m)
            tell m to ((its an) - (its gn)) < tolerance
        end |λ|
    end script
 
    script nextRefinement
        on |λ|(m)
            tell m
                set {an, gn} to {its an, its gn}
                {an:(an + gn) / 2, gn:(an * gn) ^ 0.5}
            end tell
        end |λ|
    end script
 
    an of |until|(withinTolerance, ¬
        nextRefinement, {an:(a + g) / 2, gn:(a * g) ^ 0.5})
end agm
 
-- TEST ----------------------------------------------------------------------
on run
 
    agm(1, 1 / (2 ^ 0.5))
 
    --> 0.847213084835
 
end run
 
-- GENERIC FUNCTIONS ---------------------------------------------------------
 
-- until :: (a -> Bool) -> (a -> a) -> a -> a
on |until|(p, f, x)
    set mp to mReturn(p)
    set v to x
    tell mReturn(f)
        repeat until mp's |λ|(v)
            set v to |λ|(v)
        end repeat
    end tell
    return v
end |until|
 
-- Lift 2nd class handler function into 1st class script wrapper 
-- mReturn :: Handler -> Script
on mReturn(f)
    if class of f is script then
        f
    else
        script
            property |λ| : f
        end script
    end if
end mReturn