A naturally occurring prime generator

Question

There are quite a large number of prime generating functions. Pretty much all of them are constructed and are based on the sieve of Eratosthenes, the Möbius function or the Wilson's theorem and are generally infeasible to compute in practice. But there are also generators, that have a very easy structure and were found by accident.

In 2003 Stephen Wolfram explored a class of nested recurrence equations in a live computer experiment at the NKS Summer School. A group of people around Matthew Frank followed up with additional experiments and discovered an interesting property of the simply recurrence

a(n) = a(n-1) + gcd(n,a(n-1))

with the start value of a(1) = 7. The difference a(n) - a(n-1) = gcd(n,a(n-1)) always seemed to be 1 or a prime. The first few differences are (OEIS A132199):

1, 1, 1, 5, 3, 1, 1, 1, 1, 11, 3, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 23, 3, 1, 1, 
1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 47, 3, 1, 5, 3, ...

If we only omit the 1s we get the following sequence (OEIS A137613):

5, 3, 11, 3, 23, 3, 47, 3, 5, 3, 101, 3, 7, 11, 3, 13, 233, 3, 467, 3, 5, 3, 
941, 3, 7, 1889, 3, 3779, 3, 7559, 3, 13, 15131, 3, 53, 3, 7, 30323, 3, ...

Eric S. Rowland proved the primeness of each element in this list a few years later. As you can see, the primes are mixed and some of them appear multiple times. It also has been proven, that the sequence includes infinitely many distinct primes. Furthermore it is conjectured, that all odd primes appear.

Because this prime generator was not constructed but simply found by accident, the prime generator is called "naturally occurring". But notice that in practice this generator is also quite infeasible to compute. As it turns out, a prime p appears only after (p–3)/2 consecutive 1s. Nevertheless implementing this prime generator will be your task.

Challenge:

Write a function or a program that prints the first n elements of the sequence A137613 (the sequence without the 1s). You can read the input number n >= 0 via STDIN, command-line argument, prompt or function argument. Output the first n elements in any readable format to STDOUT, or return an array or a list with these values.

This is code-golf. Therefore the shortest code wins.

Leaderboard:

Here is a Stack Snippet to generate both a regular leaderboard and an overview of winners by language. To make sure that your answer shows up, please start your 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

var QUESTION_ID=55272;function answersUrl(e){return"http://api.stackexchange.com/2.2/questions/"+QUESTION_ID+"/answers?page="+e+"&pagesize=100&order=desc&sort=creation&site=codegolf&filter="+ANSWER_FILTER}function getAnswers(){jQuery.ajax({url:answersUrl(page++),method:"get",dataType:"jsonp",crossDomain:!0,success:function(e){answers.push.apply(answers,e.items),e.has_more?getAnswers():process()}})}function shouldHaveHeading(e){var a=!1,r=e.body_markdown.split("\n");try{a|=/^#/.test(e.body_markdown),a|=["-","="].indexOf(r[1][0])>-1,a&=LANGUAGE_REG.test(e.body_markdown)}catch(n){}return a}function shouldHaveScore(e){var a=!1;try{a|=SIZE_REG.test(e.body_markdown.split("\n")[0])}catch(r){}return a}function getAuthorName(e){return e.owner.display_name}function process(){answers=answers.filter(shouldHaveScore).filter(shouldHaveHeading),answers.sort(function(e,a){var r=+(e.body_markdown.split("\n")[0].match(SIZE_REG)||[1/0])[0],n=+(a.body_markdown.split("\n")[0].match(SIZE_REG)||[1/0])[0];return r-n});var e={},a=1,r=null,n=1;answers.forEach(function(s){var t=s.body_markdown.split("\n")[0],o=jQuery("#answer-template").html(),l=(t.match(NUMBER_REG)[0],(t.match(SIZE_REG)||[0])[0]),c=t.match(LANGUAGE_REG)[1],i=getAuthorName(s);l!=r&&(n=a),r=l,++a,o=o.replace("{{PLACE}}",n+".").replace("{{NAME}}",i).replace("{{LANGUAGE}}",c).replace("{{SIZE}}",l).replace("{{LINK}}",s.share_link),o=jQuery(o),jQuery("#answers").append(o),e[c]=e[c]||{lang:c,user:i,size:l,link:s.share_link}});var s=[];for(var t in e)e.hasOwnProperty(t)&&s.push(e[t]);s.sort(function(e,a){return e.lang>a.lang?1:e.lang<a.lang?-1:0});for(var o=0;o<s.length;++o){var l=jQuery("#language-template").html(),t=s[o];l=l.replace("{{LANGUAGE}}",t.lang).replace("{{NAME}}",t.user).replace("{{SIZE}}",t.size).replace("{{LINK}}",t.link),l=jQuery(l),jQuery("#languages").append(l)}}var ANSWER_FILTER="!t)IWYnsLAZle2tQ3KqrVveCRJfxcRLe",answers=[],page=1;getAnswers();var SIZE_REG=/\d+(?=[^\d&]*(?:&lt;(?:s&gt;[^&]*&lt;\/s&gt;|[^&]+&gt;)[^\d&]*)*$)/,NUMBER_REG=/\d+/,LANGUAGE_REG=/^#*\s*([^,]+)/;

body{text-align:left!important}#answer-list,#language-list{padding:10px;width:290px;float:left}table thead{font-weight:700}table td{padding:5px}

<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="answer-list"> <h2>Leaderboard</h2> <table class="answer-list"> <thead> <tr><td></td><td>Author</td><td>Language</td><td>Size</td></tr></thead> <tbody id="answers"> </tbody> </table></div><div id="language-list"> <h2>Winners by Language</h2> <table class="language-list"> <thead> <tr><td>Language</td><td>User</td><td>Score</td></tr></thead> <tbody id="languages"> </tbody> </table></div><table style="display: none"> <tbody id="answer-template"> <tr><td>{{PLACE}}</td><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>

Answer 1

Pyth, 14 13 bytes

meaYhP+sY-5dQ

Uses a(n) = Lpf(6 - n + sum(a(i) for i in range(1, n)) where Lpf means least prime factor.

Try it here online.

Answer 2

Python 3.5.0b1+, 95 93 bytes

Link to Python 3.5.0b1+ release

import math
def f(k,n=2,a=7,L=[]):x=math.gcd(n,a);return k and f(k-1%x,n+1,a+x,L+1%x*[x])or L

A direct implementation of the recurrence, featuring:

• Our good friend 1%x, and
• math.gcd, as opposed to fractions.gcd.

Answer 3

Julia, 110 bytes

n->(a(n)=(n≥1&&(n==1?7:a(n-1)+gcd(n,a(n-1))));i=2;j=0;while j<n x=a(i)-a(i-1);x>1&&(j+=1;println(x));i+=1end)

Ungolfed:

function a(n::Int)
    n ≥ 1 && (n == 1 ? 7 : a(n-1) + gcd(n, a(n-1)))
end

function f(n::Int)
    i = 2;
    j = 0;
    while j < n
        x = a(i) - a(i-1)
        if x > 1
            j += 1
            println(x)
        end
        i += 1
    end
end

Answer 4

PHP, 101 96 99 98 77 72 bytes

<?for(;2>$t=gmp_strval(gmp_gcd(~++$i,7+$e+=$t))or$argv[1]-=print"$t ";);

Usage:
Call the Script with an argument: php -d error_reporting=0 script.php 30
If you want to test it you need to uncomment ;extension=php_gmp.dll in your php.ini
--> extension=php_gmp.dll
Should I add the extension to my byte count? Any thoughts?

Log:
Saved 3 bytes thanks to Ismael Miguel.
Saved 26 bytes thanks to primo.

Answer 5

Pyth, 30 bytes

Very badly golfed, can be considerably reduced. Defines recursive function at front, filters .first-n, and then maps the difference.

L?tb+KytbibK7m-yhdyd.ft-yhZyZQ

Try it here online.

Answer 6

Julia, 69 67 bytes

n->(i=1;a=7;while n>0 x=gcd(i+=1,a);a+=x;x>1&&(n-=1;println(x))end)

This is a simple iterative solution to the problem. x is the difference (which is the gcd), and then I update a by adding x.

Answer 7

Haskell, 51 bytes

d=zipWith gcd[2..]$scanl(+)7d
f=($filter(>1)d).take

Note that f is a function that will return the first n elements.

Rather than computing the a(n) and then working out the differences, we compute the differences d(n) and sum them together to get a(n). (Those unfamiliar with Haskell may protest that we need a(n) first in order to get d(n), but of course lazy evaluation gets us around this problem!)

Ungolfed:

a = scanl (+) 7 d        -- yielding a(n) = 7 + d(1) + d(2) + ... + d(n-1)
d = zipWith gcd [2..] a  -- yielding d(n) = gcd(n+1, a(n))

f n = take n $ filter (> 1) d -- get rid of 1s and take the first n

Answer 8

JavaScript (ES6), 91

Recursive gcd, iterative main function. Not so fast.

Usual note: test running the snippet on any EcmaScript 6 compliant browser (notably not Chrome not MSIE. I tested on Firefox, Safari 9 could go)

F=m=>{
  for(G=(a,b)=>b?G(b,a%b):a,o=[],p=7,n=1;m;d>1&&(o.push(d),--m))
    p+=d=G(++n,p);
  return o
}

O.innerHTML=F(+I.value)

<input id=I value=10><button onclick='O.innerHTML=F(+I.value)'>-></button>
<pre id=O></pre>

Answer 9

Haskell, 74 71 66 bytes

f=($filter(>1)$tail>>=zipWith(-)$scanl(\x->(x+).gcd x)7[2..]).take

Used the trick here: https://codegolf.stackexchange.com/a/39730/43318, and made point-free.

(Previous: 71 bytes)

a=scanl(\x->(x+).gcd x)7[2..]
f m=take m$filter(>1)$zipWith(-)(tail a)a

First make the sequence of a's, and then take the differences.

(Previous: 74 bytes)

f m=take m$filter(>1)$map snd$scanl(\(x,d)->(\y->(x+y,y)).gcd x)(7,1)[2..]

Standard list functions, plus clever use of lambda function. Note this is 1 byte shorter than the more obvious

g m=take m$filter(>1)$map snd$scanl(\(x,d)n->(x+gcd x n,gcd x n))(7,1)[2..]

If we don't count imports, I can get this down to 66.

import Data.List
h m=take m$filter(>1)$snd$mapAccumL(\x->(\y->(x+y,y)).gcd x)7[2..]

Answer 10

PARI/GP, 60 bytes

a(n)=a=7;i=1;while(n,if(1<d=gcd(i+=1,a),n-=1;print(d));a+=d)

Taken more or less straight from the definition a(n) - a(n-1) = gcd(n, a(n-1))

Output for a(20):

5
3
11
3
23
3
47
3
5
3
101
3
7
11
3
13
233
3
467
3

Answer 11

C++, 193 182 180 172 bytes

Thanks @Jakube - saved 8 bytes on output.

int g(int a,int b){return a==b?a:a>b?g(b,a-b):g(a,b-a);}void f(int *r,int n){int m=1,i=0,a=7,b;while(i<n){b=g(a,++m);if(b>1){r[i]=b;++i;}a+=b;}}int main(){int r[6];f(r,6);}

Answer 12

Mathematica, 59 bytes

For[i=j=1;a=7,i<=#,,d=GCD[++j,a];If[d>1,Print@d;i++];a+=d]&