# This crude Maple program is based on Liskovec's (1970) paper. Integer partitions of n are written in frequency or multiplicity notation; 
# e.g., partition 1+2+2 of 5 is written as [1,2,0,0,0] because 1*1 + 2*2 + 3*0 + 4*0 + 5*0 = 5.

# This procedure calculates the function pi of an integer partition in Liskovec (1970, p. 39).

pip := proc(a) local n, w, i: n := nops(a): w := 1: for i to n do w := w*i^a[i]: end do: w: end proc:

# This procedure calculates the factorial of an integer partition (p. 39 in Liscovec).

fac := proc(a) local n, w, i: n := nops(a): w := 1: for i to n do w := a[i]!*w: end do: w: end proc:

# This procedure calculates the function K of an integer partition in Liskovec (1970, p. 42).

kap := proc(a) local n, i: n := nops(a): add(a[i], i = 1 .. n): end proc:

# This procedure calculates the sum of two integer partitions a and b that might have different lengths.

adp := proc(a, b) local na, nb, c, i: na := nops(a): nb := nops(b): if nb < na then c := a: else c := b: end if: 
       for i to min(na, nb) do c[i] := a[i] + b[i]: end do: c: end proc:

# This procedure transforms an integer partition p of n from standard notation to frequency notation (courtesy of Peter Luschny).

trans := proc(p, n) local m, j: [seq(add(`if`(p[j] = m, 1, 0), j = 1 .. nops(p)), m = 1 .. n)]: end proc:

# This procedure calculates highest exponent of 2 that divides a number n (courtesy of Peter Luschny); viz. A007814. 
# Liskovec (1970, p. 43) uses 2^A007814(n), but that is not necessary.

expo2 := proc(n) padic[ordp](n, 2): end proc:

# This procedure calculates the polynomial p0(a,x) of an integer partition a (see p. 40 in Liskovec).

p0 := proc(a, x) local n, w, s, m: n := nops(a): w := 1: 
      for s to n - 1 do for m from s + 1 to n do w := (1 + x^lcm(s, m))^(gcd(s, m)*a[s]*a[m])*w: end do: end do: 
      for s to n do w := w*(1 + x^s)^(s*binomial(a[s], 2) + floor(1/2*s - 1/2)*a[s]): end do: 
      for s to n do if 0 = s mod 2 then w := w*(1 + x^(1/2*s))^a[s]: end if: end do: w: end proc:

# This procedure calculates the polynomial p0(+/-)(a,b,x) of two integer partitions a and b which may not necessarily be of the same length. See p. 43 in Liskovec.

p0pm := proc(a, b, x) local na, nb, w, s, m: na := nops(a): nb := nops(b): w := 1: 
       for s to na do for m to nb do if expo2(m) < expo2(s) then w := w*(1 + x^lcm(s, m))^(gcd(s, m)*a[s]*b[m]): end if: end do: end do: 
       for s to na do for m to nb do if expo2(s) <= expo2(m) then w := w*(1 - x^lcm(s, m))^(gcd(s, m)*a[s]*b[m]): end if: end do: end do: w: end proc:

# This procedure calculates the polynomial p0(-)(a,x) of an integer partition a (see p. 43 in Liskovec).

p0m := proc(a, x) local n, w, s, m; n := nops(a): w := 1: 
       for s to n - 1 do for m from s + 1 to n do if expo2(s) = expo2(m) then w := w*(1 + x^lcm(s, m))^(gcd(s, m)*a[s]*a[m]): end if: end do: end do: 
       for s to n - 1 do for m from s + 1 to n do if expo2(s) <> expo2(m) then w := w*(1 - x^lcm(s, m))^(gcd(s, m)*a[s]*a[m]): end if: end do: end do: 
       for s to n do w := w*(1 + x^s)^(s*binomial(a[s], 2) + floor(1/2*s - 1/2)*a[s]): end do: 
       for s to n do if 0 = s mod 2 then w := w*(1 - x^(1/2*s))^a[s]: end if: end do: w: end proc:

# The following procedure calculates a term in the summation formula (5) on p. 42 in Liskovec (1970). 
# Given two integer partitions a and b (of not necessarily the same length), it is a polynomial in x.

term := proc(a, b, x) 2^(-kap(adp(a, b)))*p0(a, x)*p0pm(a, b, x)*p0m(b, x)/(fac(a)*fac(b)*pip(adp(a, b))): end proc:

# The following procedure calculates the o.g.f. E0(n,x) of the n-th row of the irregular triangle T(n,k) (see p. 42 in Liskovec).

E0 := proc(n, x) local w, n1, p1, p2; w := 0; 
      for n1 from 0 to n do for p1 in combinat:-partition(n1) do 
      for p2 in combinat:-partition(n - n1) do w := w + term(trans(p1, n1), trans(p2, n - n1), x): end do: end do: end do: expand(w): end proc:
 
# The following procedure gives the n-th row of  of the irregular triangle T(n,k).

roww := proc(n) local k: seq(coeff(E0(n, x), x, k), k = 0 .. 1/2*n*(n - 1)): end proc:

# The following procedure gives the first n rows of triangle T(n,k) in flattened form.

flattenn := proc(n) local n1: seq(roww(n1), n1 = 0 .. n): end proc: