program PNfiniteBeven
! Identification of leading pairwise-nonrecurrant sequences, by Thomas Bewley.
logical :: periodic=.false. ! Specify the periodic or nonperiodic case.
integer :: m=6, show=40     ! Specify the value of m to test and the # of results to show.
! Note: compile with, eg: gfortran PNfiniteBeven.f90 -o PNfiniteBeven
integer :: r,n,num=0,location(9)=0,pairs(63)=0,a(9)=0,answer(64)=0,shown=0,best(64)=0
integer*8 :: group(20000000)=0
! Set values for r (# of times each symbol is repeated) and n (# of total slots)
! [In the nonperiodic case, the following value of r is a good one to check for 3<m<6]
! [In the periodic case, the following value of r is a good one to check for 2<m<8]
! [For other ranges of m, look at the corresponding PNfiniteFastUneven code.]
if (periodic) then; r=m; else; r=m+2; end if; n=m*r; best(1:n)=m; n=49
! Now find those groupings of r locations in the n available slots in which a given
! symbol can be repeated without violating the properties of pairwise nonrecurrance.
write(*,'(A,L2,3(A,I2))') 'periodic=',periodic,', m=',m,', r=',r,', n=',n
call system('date');
write(*,'(2(A,I2))') 'Looking for pairwise nonrecurrent groupings of length n=',n,' with r=',r
do k=1,n; pairs(1:n-1)=100; location(1)=k; call FindGroups(2); end do;
write(*,'(A,I9,A,I2)') 'Found ',num,' pairwise nonrecurrent groupings with r=',r
call system('date');
write(*,'(A,I2,A)') 'Now looking at the possible ways to fit',m,' of these groupings together'
do loop=1,num; a(1)=loop;
  if (.not. BTEST(group(a(1)),n-1)) then
    call system('date'); write(6,*) 'Exhaustively checked all possible fits.'; stop
  end if; call Check(2)
end do
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
contains
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
recursive subroutine FindGroups(level)
  integer :: level,s1,s2,i,j
  logical :: success
  do i=location(level-1)+1,n-(r-level)
  	do j=1,n-1; if (pairs(j)>=level) pairs(j)=100; end do
	location(level)=i; success=.true.
	l1: do j=level-1,1,-1
	  s1=i-location(j); s2=location(j)-i+n
	  if (periodic) then
        if (pairs(s1)>level .and. pairs(s2)>level) then; pairs(s1)=level; pairs(s2)=level
        else; success=.false.; exit l1; end if
      else
        if (pairs(s1)>level) then; pairs(s1)=level
		else; success=.false.; exit l1; end if
      end if
	end do l1
	if (success) then
	  if (level==r) then;
		  num=num+1; do j=1,r; group(num)=IBSET(group(num),n-location(j)); end do
!		  write(6,*) num,location,group(num); call FLUSH(6)
	  else; call FindGroups(level+1); end if
	end if
  end do
end subroutine FindGroups
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
recursive subroutine Check(level)
logical   :: test
integer   :: pos,x,x1,x2
integer*8 :: temp,temp1
temp=group(a(1)); do j=2,level-1; temp=ior(group(a(j)),temp); end do
l4: do pos=n-1,0,-1; if (.not. btest(temp,pos)) exit l4; end do l4
temp=0; temp=IBSET(temp,pos);    ! note: pos is the most significant bit of the next group.
temp1=0; temp1=IBSET(temp1,pos+1);  
x1=a(level-1); x2=num ! Use bisection algorithm to find starting point in search for next
do while (x2-x1>1)    ! group (find first entry such that pos+1 is the most significant bit)
   x=(x2+x1)/2; if (temp1>group(x) .and. temp1 > 0) then; x2=x; else; x1=x; end if
end do                ! need to check sign above only for n=64, integer*8 case

l3: do i=x1,num; a(level)=i; test=.true.
  l2: do j=level-1,1,-1
    test=(test .and. (iand(group(a(j)),group(a(level)))==0)); if (.not. test) exit l2
  end do l2                        ! need to check sign below only for n=64, integer*8 case
  if (temp >= group(a(level)) .and. group(a(level)) > 0) exit l3  
  if (test) then
    if (level<m) then
      call Check(level+1)
    else
      do l=1,m; do k=1,n; if (BTEST(group(a(l)),k-1)) answer(n+1-k)=l; end do; end do
      write(6,'(100I2)') answer(1:n)-1; call FLUSH(6); shown=shown+1;
	  l5: do k=1,n;
	    if (best(k)>answer(k)) then
	      write(6,*) 'Above entry is lexographically first of all sequences found thus far'
		  best=answer; exit l5
		else if (best(k)<answer(k)) then
		  exit l5
		end if
	  end do l5
	  if (shown==show) then; call system('date'); stop; end if
    end if
  end if
end do l3
end subroutine check
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
end program PNfiniteBeven