Source code in the pari/gp programming language

\\ Creating Vlf - Vector of ludic numbers' flags,
\\ where the index of each flag=1 is the ludic number.
\\ 2/28/16 aev
ludic(maxn)={my(Vlf=vector(maxn,z,1),n2=maxn\2,k,j1);
for(i=2,n2, 
    if(Vlf[i], k=0; j1=i+1;
       for(j=j1,maxn, if(Vlf[j], k++); if(k==i, Vlf[j]=0; k=0))
      ); 
   );
return(Vlf);
}

{
\\ Required tests:
my(Vr,L=List(),k=0,maxn=25000);
Vr=ludic(maxn);
print("The first 25 Ludic numbers: ");
for(i=1,maxn, if(Vr[i]==1, k++; print1(i," "); if(k==25, break)));
print("");print("");
k=0;
for(i=1,999, if(Vr[i]==1, k++));
print("Ludic numbers below 1000: ",k);
print("");
k=0;
print("Ludic numbers 2000 to 2005: ");
for(i=1,maxn, if(Vr[i]==1, k++; if(k>=2000&&k<=2005, listput(L,i)); if(k>2005, break)));
for(i=1,6, print1(L[i]," "));
print(""); print("");
print("Ludic Triplets below 250: ");
for(i=1,250, if(Vr[i]&&Vr[i+2]&&Vr[i+6], print1("(",i," ",i+2," ",i+6,") ")));
}

\\ Creating Vl - Vector of ludic numbers.
\\ 2/28/16 aev
ludic2(maxn)={my(Vw=vector(maxn, x, x+1),Vl=Vec([1]),vwn=#Vw,i);
while(vwn>0, i=Vw[1]; Vl=concat(Vl,[i]);
 Vw=vector((vwn*(i-1))\i,x,Vw[(x*i+i-2)\(i-1)]); vwn=#Vw
); return(Vl);
}
{
\\ Required tests:
my(Vr,L=List(),k=0,maxn=22000,vrs,vi);
Vr=ludic2(maxn); vrs=#Vr;
print("The first 25 Ludic numbers: ");
for(i=1,25, print1(Vr[i]," "));
print("");print("");
k=0;
for(i=1,vrs, if(Vr[i]<1000, k++, break));
print("Ludic numbers below 1000: ",k);
print("");
k=0;
print("Ludic numbers 2000 to 2005: ");
for(i=2000,2005, print1(Vr[i]," "));
print("");print("");
print("Ludic Triplets below 250: ");
for(i=1,vrs, vi=Vr[i]; if(i==1,print1("(",vi," ",vi+2," ",vi+6,") "); next); if(vi+6<250,if(Vr[i+1]==vi+2&&Vr[i+2]==vi+6, print1("(",vi," ",vi+2," ",vi+6,") "))));
}