C Count all complete rulers of a given length (OEIS A103295) C and all complete rulers with a given number of segments (A103296) C C Author: Hugo Pfoertner http://www.pfoertner.org/ C URL of this file: C http://www.randomwalk.de/sequences/a103295.txt C C Version history C Mar 4 2005 Difference counting by shifts C Mar 3 2005 Initial version C C This program depends on the availability of 64-bit integers: implicit integer*8 (a-z) logical btest C Length of lookup table for bit counting parameter(tbits=21, mask = 2**tbits-1) C mmax has to be set to (range of ruler lengths to be investigated)+1 parameter (mmax=37) dimension bitset(0:mask), rulcnt(0:mmax) C Prepare lookup table for bit counting do 20 i = 0, mask bc = 0 do 25 k = 0,tbits 25 if ( btest(i,k) ) bc = bc + 1 bitset(i) = bc 20 continue C array for counting rulers do 15 i = 0, mmax rulcnt(i) = 0 15 continue C C Loop over number of interior positions do 100 m = 1, mmax last = 2**m-1 count = 0 C Loop over all possible rulers with length m+1 do 10 r = 0, last C set all bits (no difference covered) difmis = last C C Check differences of set bits C Distances of set bits from left boundary C This clears the bits at all positions of the marks difmis = iand ( difmis, not(-not(r)-1) ) C...+....1....+....2....+....3....+....4....+....5....+....6....+....7.. do 30 j = 1, m C Clear bit positions measured from right boundary if ( btest(r,m-j) ) difmis = ibclr ( difmis, m-j ) C Check only positions not cleared so far to avoid the expensive shift if (btest(difmis,j)) then C If there is at least one common position in the original bit field and C the bit field shifted by j positions, then difference j is cleared if ( iand( r, ishft(r,j)) .ne. 0 ) difmis=ibclr ( difmis, j ) endif 30 continue C check if all bits have been cleared if ( difmis .eq. 0 ) then count = count + 1 C number of set bits in this ruler bits = bitset(iand(r,mask))+ bitset(ishft(r,-tbits)) rulcnt(bits) = rulcnt(bits) + 1 endif C End of loop over all rulers of length m 10 continue write (*,*) m, count 100 continue do 200 i = 1, mmax if ( rulcnt(i) .gt. 0 ) write(*,*) i, rulcnt(i) 200 continue end