[seqfan] Re: I need a name for this sequence

[Aai]

By using bounds of a for every c >= 2 and calculating c_max for N<L, these 
numbers can be generated very fast with the following PARI/GP program. I use PL 
J (http://www.jsoftware.com/) to investigate these kind of problems. But AFAIU 
most of you can't read/understand such code, so I translated the method into a 
PARI/GP program.

Be aware: I'm not an expert PARI/GP user.


Only use for N >= 15

aseq(N) = {
  local(r,cm,v,k);
  r = real(polroots (4*x^2 - x - N));
  if(r[1]>0,cm=r[1],cm=r[2]); cm=floor(cm);
  v=vector(2^floor(log(N))); k=0;
  for(c=2,cm, for(a=c+1,2*c-1,i=a-c-1;b=c*(c+i)/(1+i); 
if(b==floor(b),t=a*b+c;if(t<N,v[k++]=t))));
  nub(vecsort(vecextract(v,concat("..",k))))
}

I had to write a nub function because I use PARI/GP 2.3.5:

nub (V) = {
  local(l,W,k);
  l= length(V); W=vector(l); k=0;
  for(i=1,l-1,if(V[i]!=V[i+1],W[k++]=V[i]));
  W[k++]=V[l];
  return (vecextract(W,concat("..",k)))
}

Example (hand delimited) output:

(19:11) gp > aseq(10000)
time = 8 ms.
%5 = [14, 33, 39, 60, 64, 84, 95, 110, 138, 150, 155, 174, .. ]



-- 
Met vriendelijke groet,
=@@i