(48) (49) (50) Moebius transforms

[John Keith]

Here are three Moebius transform related programs. They are most useful for integer sequences related to the divisor function, the Moebius Mu function, and the Euler Phi function.

First, the Moebius transform. The program takes a list of integers and returns a list of integers of the same length. It requires ListExt and MOB, Gerald H's program for Mu(n) detailed in this thread.

Code:


\<< DUP SIZE R\->I \-> n
  \<< 1 n
    FOR k k DIVIS DUP2 LPICK SWAP REVLIST Mu * LSUM SWAP
    NEXT DROP n \->LIST
  \>>
\>>

Next, the inverse Moebius transform. Also requires ListExt.

Code:


\<< DUP SIZE R\->I \-> n
  \<< 1 n
    FOR k DUP k DIVIS LPICK LSUM SWAP
    NEXT DROP n \->LIST
  \>>
\>>

Finally, a program for Dirichlet convolution, a closely related operation. Note that this is not a convolution of two lists, but rather of two functions over a list of the natural numbers.

The program requires two programs for the desired functions on levels 2 and 3, and an integer representing the length of the list to be returned on level 1.

Code:


\<< \-> f g n
  \<< 1 n
    FOR k k DIVIS DUP f EVAL SWAP REVLIST g EVAL * 0 + \GSLIST
    NEXT n \->LIST
  \>>
\>>

For example, if level 3 contains \<< MOB \>> from the second link above
and level 2 contains \<< EULER \>> and level 1 has the number 20,
the program will return the first 20 terms of A007431.

For anyone interested in this area of number theory, I also highly recommend Gerald H's SUMDIVISORS, detailed in this thread.

EDIT: Related to the post below, here is an alternate version of the Dirichlet convolution program which takes two sequences (lists) rather than functions. The list on level 2 must have at least as many terms as the list on level 1.

Code:


\<< DUP SIZE R\->I \-> n
  \<< 1 n
    FOR k DUP2 k DIVIS SWAP OVER LPICK UNROT LPICK REV * LSUM UNROT
    NEXT DROP2 n \->LIST
  \>>
\>>

Another note regarding the program here. If one desires to do "Dirichlet deconvolution", the method is very simple. With the result list on level 2 and one of the original lists on level 1, take the Dirichlet inverse of the level 1 list, then Dirichlet convolve the two lists.

To put it more formally, if h is the Dirichlet convolution of f and g, then g is the Dirichlet convolution of h and the Dirichlet inverse of f.

[John Keith]

Related to the programs above, here is a program for the Dirichlet inverse of a sequence. The Dirichlet inverse is defined in the Wikipedia link in the post above.

This program is defined in terms of a sequence rather than a function so that it can be applied to sequences that are not defined in terms of a simple function. The first term of the sequence must be 1, otherwise the result will be incorrect. This program also requires the ListExt Library.

Code:


\<< DUP SIZE R\->I \-> b n
  \<< { 1 } 2 n
    FOR k k DIVIS b OVER TAIL LPICK UNROT OVER SWAP REV TAIL LPICK ROT * LSUM NEG +
    NEXT
  \>>
\>>

For the reason stated above, I have added an alternate version of the Dirichlet convolution program for sequences rather than functions.

[John Keith]

Updating the Moebius transform program with a new version which is 23 bytes larger but about 30% faster. The improvement comes from pre-calculating the Moebius Mu function of all divisors.

Code:


\<< DUP SIZE R\->I \-> n
  \<< n LSEQ Mu SWAP 1 n
    FOR k DUP2 k DIVIS SWAP OVER LPICK UNROT REV LPICK * LSUM UNROT
    NEXT DROP2 n \->LIST
  \>>
\>>