

 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


with Ada.Text_IO, Ada.Numerics.Generic_Elementary_Functions;
 
procedure Arith_Geom_Mean is
 
   type Num is digits 18; -- the largest value gnat/gcc allows
   package N_IO is new Ada.Text_IO.Float_IO(Num);
   package Math is new Ada.Numerics.Generic_Elementary_Functions(Num);
 
   function AGM(A, G: Num) return Num is
      Old_G: Num;
      New_G: Num := G;
      New_A: Num := A;
   begin
      loop
         Old_G := New_G;
         New_G := Math.Sqrt(New_A*New_G);
         New_A := (Old_G + New_A) * 0.5;
         exit when (New_A - New_G) <= Num'Epsilon;
         -- Num'Epsilon denotes the relative error when performing arithmetic over Num
      end loop;
      return New_G;
   end AGM;
 
begin
   N_IO.Put(AGM(1.0, 1.0/Math.Sqrt(2.0)), Fore => 1, Aft => 17, Exp => 0);
end Arith_Geom_Mean;