This is a front-end for the Online Encyclopedia of Integer Sequences, made by Christian Perfect. The idea is to provide OEIS entries in non-ancient HTML, and then to think about how they're presented visually. The source code is on GitHub.
with(group); permul := (a,b) -> mulperms(b,a); swap := (p,i,j) -> convert(permul(convert(p,'disjcyc'),[[i,j]]),'permlist',nops(p)); PermRank3Aux := proc(n, p, q) if(1 = n) then RETURN(0); else RETURN((n-p[n])*((n-1)!) + PermRank3Aux(n-1,swap(p,n,q[n]),swap(q,n,p[n]))); fi; end; PermRank3R := p -> PermRank3Aux(nops(p),p,convert(invperm(convert(p,'disjcyc')),'permlist',nops(p))); PermRank3L := p -> PermRank3Aux(nops(p),convert(invperm(convert(p,'disjcyc')),'permlist',nops(p)),p); # a(n) = PermRank3L(PermUnrank3R(n)) or PermRank3R(PermUnrank3L(n)) or PermRank3L(convert(invperm(convert(PermUnrank3L(j), 'disjcyc')), 'permlist', nops(PermUnrank3L(j))))
from operator import mul from sympy import prime, factorial as f from sympy.ntheory.factor_ import core def a007623(n, p=2): return n if n0 else '0' for i in x)[::-1] return 0 if n==1 else sum([int(y[i])*f(i + 1) for i in range(len(y))]) def a(n): return 1 if n==0 else a275732(n)*a(a257684(n)) def ok(n): return 1 if n==0 else core(a(n))==a(n) print([n for n in range(201) if ok(n)]) # Indranil Ghosh, Jun 19 2017
n=14 (factorial base "210") is included because 2 occurs in position 3 and 1 occurs in position 2, thus as (3-2) = 1 <> 2, 2 does not "hit" digit 1. n=15 ("211") is NOT included because 2 occurring in position 3 hits the rightmost 1 in position 1 (as 3-2 = 1), and moreover, also the middle 1 hits the rightmost 1 as 2-1 = 1.
For n=27, which in factorial base (A007623) is "1011" and encodes (in A060118-order) permutation "23154" with one 3-cycle and one 2-cycle, the longest cycle has three elements, thus a(27) = 3.
N:= 100: # to get a(0) to a(N) M:= 0: A[0]:= 0: count:= 0: for m from 2 while count < N do P:= remove(t -> t[1]=1, combinat:-permute(m)); P:= map(t -> ListTools:-Reverse(subs([seq(i=m+1-i,i=1..m)],t)),P); R:= select(t -> max(map(nops,convert(P[t],disjcyc))) = 2, [$1..nops(P)]); for r in R do count:= count+1; A[count]:= r+M; if count = N then break fi; od: M:= M + nops(P); od: seq(A[i],i=0..count); # Robert Israel, Oct 28 2015
Comments