AGM Series Hole 1: Calculate the Arithmetic–Geometric Mean

Question

This question was inspired by this HNQ.

About the series

This question is now part of a series about the AGM method. This first post in the series will be about actually calculating the AGM. You may treat this like any other code golf challenge, and answer it without worrying about the series at all. However, there is a leaderboard across all challenges.

What is the Arithmetic–Geometric Mean

The Arithmetic–Geometric Mean of two numbers is defined as the number that repeatedly taking the arithmetic and geometric means converges to. Your task is to find this number after some n iterations.

Clarifications

• You take three numbers, a, b, n in any reasonable format.
• For n iterations, take the arithmetic and geometric mean of a and b and set those to a and b.
• For two numbers a and b, the arithmetic mean is defined as (a + b) / 2.
• The geometric mean is defined as √(a * b).
• a and b should be approaching each other.
• Then, output both a and b.
• You don't have to worry about float imprecision and such.
• This is code-golf so shortest code in bytes wins!

Test Cases

[0, [24, 6]] -> [24, 6] 
[1, [24, 6]] -> [15.0, 12.0]
[2, [24, 6]] -> [13.5, 13.416407864998739]
[5, [24, 6]] -> [13.458171481725616, 13.458171481725616]
[10, [100, 50]] -> [72.83955155234534, 72.83955155234534]
The next one is 1/Gauss's Constant:
[10, [1, 1.41421356237]] -> [1.198140234734168, 1.1981402347341683]

Leaderboard

Stolen from Martin's series.

The following snippet will generate a leaderboard across all challenges of the series.

To make sure that your answers show up, please start every answer with a headline, using the following Markdown template:

# Language Name, N bytes

where N is the size of your submission. If you improve your score, you can keep old scores in the headline, by striking them through. For instance:

# Ruby, <s>104</s> <s>101</s> 96 bytes

/* Configuration */
var QUESTION_IDs = [66068]; // Obtain this from the url
// It will be like http://XYZ.stackexchange.com/questions/QUESTION_ID/... on any question page
var ANSWER_FILTER = "!.FjwQBrX2KXuFkv6p2lChi_RjzM19";
/* App */
var answers = [], page = 1, currentQ = -1;
function answersUrl(index) {
 return "http://api.stackexchange.com/2.2/questions/" + QUESTION_IDs.join(";") + "/answers?page=" + index + "&pagesize=100&order=desc&sort=creation&site=codegolf&filter=" + ANSWER_FILTER;
}
function getAnswers() {
 $.ajax({
 url: answersUrl(page++),
 method: "get",
 dataType: "jsonp",
 crossDomain: true,
 success: function (data) {
 answers.push.apply(answers, data.items);
 if (data.has_more) getAnswers();
 else process();
 }
 });
}
getAnswers();
var SIZE_REG = /\d+(?=[^\d&]*(?:<(?:s>((?!>).)*<\/s>|((?!>).)+>)[^\d&]*)*$)/;
var NUMBER_REG = /\d+/;
var LANGUAGE_REG = /^#*\s*([^\n,]+)(?=,)/;//
function shouldHaveHeading(a) {
 var pass = false;
 var lines = a.body_markdown.split("\n");
 try {
 pass |= /^#/.test(a.body_markdown);
 pass |= ["-", "="]
 .indexOf(lines[1][0]) > -1;
 pass &= LANGUAGE_REG.test(a.body_markdown);
 } catch (ex) {}
 return pass;
}
function shouldHaveScore(a) {
 var pass = false;
 try {
 pass |= SIZE_REG.test(a.body_markdown.split("\n")[0]);
 } catch (ex) {}
 if (!pass) console.log(a);
 return pass;
}
function getAuthorName(a) {
 return a.owner.display_name;
}
function getAuthorId(a) {
 return a.owner.user_id;
}
function process() {
 answers = answers.filter(shouldHaveScore)
 .filter(shouldHaveHeading);
 answers.sort(function (a, b) {
 var aB = +(a.body_markdown.split("\n")[0].match(SIZE_REG) || [Infinity])[0],
 bB = +(b.body_markdown.split("\n")[0].match(SIZE_REG) || [Infinity])[0];
 return aB - bB
 });
 var users = {};
 answers.forEach(function (a) {
 var headline = a.body_markdown.split("\n")[0];
 var question = QUESTION_IDs.indexOf(a.question_id);
 var size = parseInt((headline.match(SIZE_REG)||[0])[0]);
 var language = headline.match(LANGUAGE_REG)[1];
 var user = getAuthorName(a);
 var userId = getAuthorId(a);
 if (!users[userId]) users[userId] = {name: user, nAnswer: 0, answers: []};
 if (!users[userId].answers[question]) {
 users[userId].answers[question] = {size: Infinity};
 users[userId].nAnswer++;
 }
 if (users[userId].answers[question].size > size) {
 users[userId].answers[question] = {size: size, link: a.share_link}
 }
 });
 
 
 var sortedUsers = [];
 for (var userId in users)
 if (users.hasOwnProperty(userId)) {
 var user = users[userId];
 user.score = 0;
 user.completedAll = true;
 for (var i = 0; i < QUESTION_IDs.length; ++i) {
 if (user.answers[i])
 user.score += user.answers[i].size;
 else
 user.completedAll = false;
 }
 sortedUsers.push(user);
 } 
 
 sortedUsers.sort(function (a, b) {
 if (a.nAnswer > b.nAnswer) return -1;
 if (b.nAnswer > a.nAnswer) return 1;
 return a.score - b.score;
 });
 
 var place = 1;
 for (var i = 0; i < sortedUsers.length; ++i) {
 var user = sortedUsers[i];
 var row = '<tr><td>'+ place++ +'.</td><td>'+user.name+'</td>';
 for (var j = 0; j < QUESTION_IDs.length; ++j) {
 var answer = user.answers[j];
 if (answer)
 row += '<td><a href="'+answer.link+'">'+answer.size+'</a></td>';
 else
 row += '<td class="missing"></td>';
 }
 row += '<td></td>';
 if (user.completedAll)
 row += '<td class="total">'+user.score+'</td>';
 else
 row += '<td class="total missing">'+user.score+'</td>';
 row += '</tr>';
 $("#users").append(row);
 }
}

body { text-align: left !important}
#leaderboard {
 width: 500px; 
}
#answer-list {
 padding: 10px;
 width: 290px;
 float: left;
}
#language-list {
 padding: 10px;
 width: 290px;
 float: left;
}
table thead {
 font-weight: bold;
}
table td {
 padding: 5px;
}
td.total {
 font-weight: bold;
 text-align: right;
}
td.missing {
 background: #bbbbbb;
}

<script src="https://ajax.googleapis.com/ajax/libs/jquery/2.1.1/jquery.min.js"></script>
<link rel="stylesheet" type="text/css" href="//cdn.sstatic.net/codegolf/all.css?v=83c949450c8b">
<div id="leaderboard">
 <h2>Leaderboard</h2>
 <p>
 Missing scores are shown as grey cells. A grey total indicates that the user has not participated in all challenges and is not eligible for the overall victory yet.
 </p>
 <table class="_user-list">
 <thead>
 <tr><td></td><td>User</td>
 <td><a href="http://codegolf.stackexchange.com/q/66068">#1</a></td>
 <td></td><td>Total</td>
 </tr>
 </thead>
 <tbody id="users">
 </tbody>
 </table>
</div>
<table style="display: none">
 <tbody id="answer-template">
 <tr><td>{{NAME}}</td><td>{{LANGUAGE}}</td><td>{{SIZE}}</td><td><a href="{{LINK}}">Link</a></td></tr>
 </tbody>
</table>
<table style="display: none">
 <tbody id="language-template">
 <tr><td>{{LANGUAGE}}</td><td>{{NAME}}</td><td>{{SIZE}}</td><td><a href="{{LINK}}">Link</a></td></tr>
 </tbody>
</table>

• xnor: Are the initial numbers positive integers?
• Doorknob: "both a or b"—well, which one? Both, or either one?
• Maltysen: @Doorknob -_- Its both.
• Maltysen: @xnor no. Look at the last test-case.
• Reto Koradi: Making this part of a series causes a kind of unfortunate situation. This is so simple that solutions are all going to look quite similar. And posting similar solutions in languages that were already used is generally frowned upon. I wrote my solution in about 2 minutes, but it's in a language that was already used, and it's the same length. If I follow typical posting etiquette, I won't be able to participate in the series.

Answer 1

CJam, 16 bytes

{{_2$+2/@@*mq}*}

Takes input on the stack as a b n where a and b are doubles. Online demo

Answer 2

Dyalog APL, (削除) 22 (削除ここまで) (削除) 21 (削除ここまで) 15 bytes

.5∘(×ばつ.*⍨)⍣⎕

Takes (a,b) as right argument, and prompts for n:

(

×ばつ dot product of 0.5 and the right argument

, followed by

×ばつ.*⍨ "dot power" of the right argument and 0.5^*

)⍣⎕ applied numeric-prompt times.

^* "dot power" is like dot product, but using multiplication and power instead of plus and multiplication, as follows:

_n
A ×ばつ.*⍨ B is ∏ B_i^A = ∏ B₁^AB₂^A
i=1

-3 bytes thanks to ngn.

---

Old version:

{((+/÷≢)×ばつ/)⍣⍺⊢⍵}

Takes n as left argument and a b as right argument.

⊢⍵ On the RightArg
(...)⍣⍺ recalculate LeftArg times
(+/÷≢) sum divided by tally
, followed by
×ばつ/ the square root of the product.

All the test cases:

 f←{((×ばつ+/)×ばつ/)⍣⍺⊢⍵}
 0 1 2 5 10 10 f ̈ (24 6)(24 6)(24 6)(24 6)(100 50)(1,2*.5)
┌────┬─────┬────────────────┬───────────────────────┬───────────────────────┬───────────────────────┐
│24 6│15 12│13.5 13.41640786│13.45817148 13.45817148│72.83955155 72.83955155│1.198140235 1.198140235│
└────┴─────┴────────────────┴───────────────────────┴───────────────────────┴───────────────────────┘

• lirtosiast: Is f⍣⍺⊢⍵ or similar an idiom you use professionally?
• Adám: @ThomasKwa Yes, see e.g. Of⍣core⊢TREE at miserver.dyalog.com (click the big "D" and scroll to line [266]).

Answer 3

TI-BASIC, 22 bytes

Input N
For(I,1,N
{mean(Ans),√(prod(Ans
End
Ans

Does exactly what the algorithm says. Takes N from the prompt, and A and B through Ans as a two-element list.

If N is 0, the For( loop is skipped entirely.

Answer 4

JavaScript ES7, (削除) 48 (削除ここまで) 43 bytes

-5 thanks to Downgoat!

f=(n,a,b)=>n?f(n-1,(a+b)/2,(a*b)**.5):[a,b]

Very simple recursive function.

• Downgoat: (a*b)**.5 is shorter than Math.sqrt(a*b). example
• Conor O'Brien: @Downgoat That's ES7, but meh.

Answer 5

MATLAB/Octave, (削除) 69 (削除ここまで) 65 bytes

function [a,b]=r(a,b,n)
for i=1:n;j=(a+b)/2;b=(a*b)^.5;a=j;end

• Ander Biguri: You could do b=(a*b).^5 directly as you are not reusing b again in that iteration and save 4 bytes.

Answer 6

Jelly, 9 bytes

SH;P1⁄2\ðṛ¡

Try it online!

How it works

SH;P1⁄2\ðṛ¡ Input: x (vector) -- y (repetitions)
SH Take the sum (S) of x and halve (H) the result.
 P1⁄2 Take the product (P) of x and the square root (1⁄2) of the result.
 \ Combine the last two instructions in a dyadic chain.
 ; Concatenate the results to the left and to the right.
 ð Push the preceding, variadic chain; begin a new, dyadic chain.
 ṛ Return the right argument (y).
 ¡ Repeat the pushed chain y times.

Answer 7

C++, (削除) 108 (削除ここまで) (削除) 102 (削除ここまで) 100 bytes

Thank you to @RetoKoradi and @AlexA for saving me 6 bytes.

This is non-competitive, because C++ is not a good golfing language. Did this for fun :)

#include<cmath>
std::string f(float a,float b,int n){return n==0?a+" "+b:f((a+b)/2,sqrt(a*b),n-1);}

This is a simple recursion function, very similar to the JS answer.

• Reto Koradi: You can get rid of the spaces after the commas. Also, using float instead of double is shorter.
• Alex A.: You can also remove the space in the #include line.
• TheCoffeeCup: Wow, I'm stupid not to notice that. Thanks!
• nwp: I would consider f(float*s) which takes a pointer to 3 floats to be in a "reasonable format". Not sure if that actually makes it shorter.

Answer 8

K5, 15 bytes

Very literal:

{(+/x%2;%*/x)}/

In action:

 {(+/x%2;%*/x)}/[0; 24 6]
24 6
 {(+/x%2;%*/x)}/[5; 24 6]
1.345817e1 1.345817e1

Unfortunately, this does not work in oK because that interpreter does not currently support projection (currying) of adverbs. Works in the real k5.

In oK, it would currently be necessary to wrap the definition in a lambda:

 {x{(+/x%2;%*/x)}/y}[5; 24 6]
13.4582 13.4582

Answer 9

J, (削除) 18 (削除ここまで) 13 bytes

-:@+/,%:@*/^:

Usage:

 agm =: -:@+/,%:@*/^:
 5 agm 24 6
13.4582 13.4582

• randomra: Wow, this works. Conjunctions are weird. I would expect this expression to be an adverb (which it can be), but if presented with arguments it is also a function.

Answer 10

Seriously, 11 bytes

,p`;π√@æk`n

Hex Dump:

2c70603be3fb40916b606e

Try it online

Explanation:

, Read in the list as [n,a,b]
 p pop list to yield: n [a,b]
 ` `n Push a quoted function and run it n times.
 ; Duplicate [a,b] pair
 π√ Compute its product and square root it (GM)
 @ Swap the other copy of the pair to the top
 æ Compute its mean.
 k Compile the stack back into a list.

Answer 11

Japt, 24 bytes (削除) 25 33 (削除ここまで)

Saved 9 (削除) 7 (削除ここまで) bytes thank to @ETHproductions

Uo r@[VW]=[V+W /2(V*W q]

Takes advantage of ES6 destructuring.

Try it online

Ungolfed && Explanation

Uo r@[VW]=[V+W /2(V*W q]
 // Implicit: U: 1st input, V: 2nd input, W: 3rd input
Uo // Range from 0 to 1st input
r@ // Loop over range
 [V,W]= // Set 2nd and 3rd input to...
 [V+W /2, // Add 2nd and 3rd inputs, divide by 2
 (V*W q] // Multiple 2nd and 3rd inputs, find square root
 // Set's to the above respectively 
 // Implicit: return [V,W]

• ETHproductions: Uo generates a range of numbers from 0 to U, so Uo m@[V,W]=[V+W /2,(V*W q] should work. (Untested)
• ETHproductions: Oh, and you shouldn't need the commas at all. :)
• ETHproductions: Oh dear, this fails for any U other than 1, outputting each loop as it goes. Here's one that works properly: Uo £[VW]=[V+W /2(V*W q]};[VW]
• Downgoat: @ETHproductions thanks, but using r seemed to work also
• ETHproductions: Ah, r works much better; nice job outsmarting me in my own language. ;)

Answer 12

Matlab, 54 bytes

function x=f(x,n)
for k=1:n
x=[mean(x) prod(x)^.5];end

Example:

>> f([24 6], 2)
ans =
 13.500000000000000 13.416407864998739

Answer 13

Pyth, 12

u,.OG@*FG2EQ

Test Suite

Explanation

u,.OG@*FG2EQ ## implicit: Q = eval(input())
u EQ ## reduce eval(input()) times, starting with Q
 ## the reduce lambda has G as the previous value and H as the next
 .OG ## arithmetic mean of last pair
 @*FG2 ## geometric mean of last pair, uses *F to get the product of the list
 ## and @...2 to get the square root of that
 , ## join the two means into a two element list

• orlp: Forgot about @ and .O, but I didn't even know the new purpose of E.
• FryAmTheEggman: @orlp ah, didn't see your post, my bad I would have just suggested these in the comments. And yeah, keeping track of all the changing stuff is a bit of a struggle :P

Answer 14

Minkolang v0.14, 23 bytes

Try it here!

$n[$d+2$:r*1Mi2%?!r]$N.
$n C get all input C
 [ ] C pop N; repeat inner N times C
 $d C duplicate stack [1,2] => [1,2,1,2] C
 + C add top two elements C
 2$: C divide by two C
 r C reverse stack (get the other two) C
 * C multiply them together C
 1M C take square root C
 i2%?!r C reverse the stack if an odd step number C
 $N C output stack
 1M C take square root C
 i C get step in for loop C

Answer 15

Pyth, 15 bytes

u,^*FG.5csG2vzQ

Answer 16

Python 3, (削除) 65 (削除ここまで) 55 bytes

Thanks to mathmandan for pointing out a shorter version using the lambda operator.

f=lambda a,b,n:f((a+b)/2,(a*b)**.5,n-1)if n else(a,b)

My original version:

def f(a,b,n):
 if n:f((a+b)/2,(a*b)**.5,n-1)
 else:print(a,b)

To my chagrin, a recursive function (a la the JavaScript and C++ answers) was shorter than a simple for loop.

• mathmandan: You can shorten this a little with lambda and the ternary if/else operator: f=lambda a,b,n:f((a+b)/2,(a*b)**.5,n-1)if n else(a,b)
• mathmandan: No problem! (Also, I think this is 53 bytes.)
• Jack Brounstein: The .py file I have saved is listed as being 55 bytes. Is there a better way to calculate program size?
• mathmandan: Sometimes people on this site copy and paste their code into mothereff.in/byte-counter . If you're wondering about the discrepancy, I would guess that Windows is inserting an unnecessary newline character at the end of your .py file (and Windows counts a newline as 2 bytes instead of 1). Either way, you don't have to count that last newline as part of your code for scoring purposes. If you do post a multiple-line entry, you should count 1 for each newline character, not 2, and not including any newline at the end of your last line of code. (As far as I understand the rules anyway!)

Answer 17

R, 66 bytes

f=function(a,b,n){while(n){x=(a+b)/2;b=(a*b)^.5;n=n-1;a=x};c(a,b)}

Usage:

> f(24,6,0)
[1] 24 6
> f(24,6,1)
[1] 15 12
> f(24,6,2)
[1] 13.50000 13.41641
> f(24,6,3)
[1] 13.45820 13.45814
> f(24,6,4)
[1] 13.45817 13.45817
> f(100,50,10)
[1] 72.83955 72.83955
> f(1,1.41421356237,10)
[1] 1.19814 1.19814

• Alex A.: You can remove the function name to save 2 bytes.

Answer 18

Mathematica, (削除) 31 (削除ここまで) 30 bytes

Saved one byte thanks to Martin Büttner.

{+##/2,(1##)^.5}&@@#&~Nest~##&

Usage:

In[1]:= {+##/2,(1##)^.5}&@@#&~Nest~##&[{24, 6}, 5]
Out[1]= {13.4582, 13.4582}

Answer 19

Lua, 62 bytes

n,a,b=...for i=1,n do a,b=(a+b)/2,math.sqrt(a*b)end print(a,b)

Uses command line arguments from ... to assign to n, a and b, a nifty trick I learned about Lua recently.

Answer 20

Haskell, 40 bytes

(!!).iterate(\(a,b)->((a+b)/2,sqrt$a*b))

An anonymous function. Example usage:

>> let f=(!!).iterate(\(a,b)->((a+b)/2,sqrt$a*b)) in f (1.0,1.41421356237) 10
(1.198140234734168,1.1981402347341683)

The lambda function (\(a,b)->((a+b)/2,sqrt$a*b)) takes the arithmetic and geometric mean on a tuple. This is iterated starting with the first input (a tuple), and then (!!) indexes the second input to specify the number of iterations.

Answer 21

Perl, 60 bytes

perl -ape'F=($F[0]/2+$F[1]/2,sqrt$F[0]*$F[1])for 1..shift@F;$_="@F"'

N.B.: Per this meta post, I believe I've got the scoring correct. The actual code (between single quotes) is 58 characters, then I added +2 for a and p flags as that's the difference from the shortest invocation, perl -e'...'

Vague complaints

I have this nagging feeling I'm missing an obvious improvement. I know, "welcome to code golf", but I mean more than usual I believe there's an easy opportunity to shorten this.

Early on, I had messed around with using $\ as the second term with some success, but the above approach ended up being 2 bytes shorter, even with the extra ap flags required. Similarly, avoiding the explicit $_ assignment would be nice, but the loop makes that difficult.

The shift@F bugs me, too; if I don't do it that way, though (or use @F=(0,...,...) instead, which doesn't save any bytes), there's an off-by-one error with the @F assignment.

Example

echo 5 24 6 | perl -ape'F=($F[0]/2+$F[1]/2,sqrt$F[0]*$F[1])for 1..shift@F;$_="@F"'

Outputs

13.4581714817256 13.4581714817256

Answer 22

Julia, 49 bytes

(a,b,n)->(for i=1:n;a,b=(a+b)/2,√(a*b)end;(a,b))

Pretty direct iterative algorithm. Using the√ symbol and the multiple return saves a few bytes, but the for loop syntax costs a few.

Answer 23

Haskell, 47 Bytes

f a b 0=(a,b)
f a b n=f((a+b)/2)(sqrt$a*b)(n-1)

• Damien: you could save some bytes taking a b as a couple in f: f x 0=x;f(a,b)n=f((a+b)/2,sqrt$a*b)$n-1
• xnor: And define the function infix.

Answer 24

Julia, 42 bytes

f(a,b,n)=n>0?f((a+b)/2,(a*b)^.5,n-1):(a,b)

This is a recursive function f that accepts three numbers and returns a tuple.

Ungolfed:

function f(a::Real, b::Real, n::Integer)
 if n > 0
 # Recurse on the arithmetic and geometric means, decrementing n
 return f((a + b) / 2, sqrt(a * b), n - 1)
 else
 # Return the pair
 return (a, b)
 end
end

Answer 25

LabVIEW, 21 LabVIEW Primitives

Primitives counted as per this meta post.

pretty staightforward not much to explain.

Answer 26

Python 2, (削除) 62 (削除ここまで) (削除) 61 (削除ここまで) 62 bytes

def f(a,b,n):
 while n:a,b=(a+b)/2.,(a*b)**.5;n-=1
 print a,b

• lirtosiast: The program should only print once, when it terminates.
• wflynny: My misunderstanding. Fixed.

Answer 27

CJam, 16 bytes

{{_:+2/\:*mq]}*}

This is an anonymous function. The input is a list with the two values (as doubles), followed by the iteration count. Try it online with I/O code for testing.

I wouldn't normally have posted this because @PeterTaylor posted an equally long CJam answer before I saw the question. But since this is advertised as the start of a series, I wanted to keep my options open in case the series is interesting.

While the length is the same as Peter's answer, the code is not. I chose a different input format by taking the two values in a list, where Peter used separate values. So while there's not much to it with either input format, the code looks quite different.

{ Start loop over number of iterations.
 _ Copy the current pair of values.
 :+ Reduce pair with + operator.
 2/ Divide by 2.
 \ Swap second copy of pair to top.
 :* Reduce pair with * operator.
 mq Calculate square root.
 ] Wrap the two new values in a list for next iteration.
}* End iteration loop.

Answer 28

><> (Fish), 66 bytes

Uses newtons method to approximate the square root in 15 steps

iii>:?v~naon;
*&f:?v1円-@:&$:@$:@+2,@
$@:$@//,2+,@:
4{~$~/\@@1-21.0

Try it

enter image description here

Answer 29

Scala, 58 bytes

Golfed version. Try it online!

(n,a,b)=>if(n>0)f(n-1,(a+b)/2,Math.sqrt(a*b))else s"$a,$b"

Answer 30

Thunno 2, 6 bytes

{[m;pƭ

Try it online!

Takes input on the stack. This default was +11/-5 at the time of posting (Wayback Machine).

Explanation

{[m;pƭ # Implicit input
{ # First input no. times:
 [ ; # Parallelly apply:
 m # Arithmetic mean
 pƭ # Square root of product
 # Implicit output